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SIX-DEGREE-OF-FREEDOM EQUATIONS OF MOTION FOR 
A MANEUVERING RE-ENTRY VEHICLE 


ABSTRACT 

Equations of motion for a re-entry vehicle are developed for use with a 
digitally operated, fixed-base simulator. The equations provide: 

1. Six-degrees-of-freedom for a vehicle of non-varying mass whose 
trajectory and orientation are given with respect to a spherical, 
rotating earth. 

2. An inverse-square gravitational force, a rotating atmosphere, and 
relevant atmospheric properties. 

3. Aerodynamic forces and moments in coefficient form, and reaction 
control moments. 

The only restriction on attitude, maneuver, or flight path is that flights directly 
over the North and South poles are excluded. 

The equations, as described, have been programmed on an IBM 704 
and specific solutions were obtained to show effects of maneuvers, pilot's 
inputs, different geophysical models, and mathematical simplifications. 

The report contains descriptive chapters on geophysics and aerodynamics 
and control which are especially pertinent to re-entry simulation. In addition, 
Appendices extend the equations for varying mass, large thrust, and earth oblate¬ 
ness with related effects. 


Reproduction of this publication 
in whole or in part is permitted 
for any purpose of the United 
States Government 
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FOREWORD 


The application of high speed general purpose type digital computers 
to real time simulation vas pioneered by the Naval Training Device Center 
and continued with the development of UDOFTT (Universal Digital Operational 
Flight Trainer Tool) under joint sponsorship with the Air Force, Wright Air 
Development Division. UDOFTT simulated a subsonic F9F and a supersonic 
F100 aircraft. The present study was initiated as a continuation of this 
development, to further explore the capabilities and requirements of real¬ 
time digital simulation in other areas. 

The scope of the study was limited to an Investigation of the equations 
of motion required for manned hypersonic flight and the possible simplifica¬ 
tion of these equations for computation in a digital computer for real-time 
simulation. A hypothetical vehicle was generated to be studied and to further 
limit the scope of the study, only the re-entry portion of the flight was 
to be considered. The study did have the advantage that simulation for 
training may be limited in frequency and response to those equations to which 
a pilot can be expected to react. 

The present study concludes that considerable simplification can be 
effected in the general sense and additional simplification may be tolerable 
in specific training areas. However, the full flight dynamics of piloted 
hypersonic or near-orbital vehicles over the complete flight regime from 
boost, through orbit, re-entry and landing are too complex for simulation 
on UDOFTT or present state-of-the-art computers. The study suggests that 
a breakdown into simpler regimes of limited complexity may be necessary to 
train for portions of the total flight. The possibility exists that high 
speed magnetic tape storage of portions of the program with entry into the 
computer during transitions from one phase to another would permit a continuous 
simulation as perceived by the pilot. 


Head, Computer Branch 

U.S. Naval Training Device Center 
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SYMBOLS AND DEFINITIONS 

Symbols, definitions, and conventions are illustrated in Figures 3. 1 
through 3. 5 and elsewhere in the report. 

AXES 


The following axes are right-handed orthogonal sets with origin at the 
center of gravity (c. g. ), and a longitudinal plane of mirror symmetry for the 
aircraft is assumed. 


Wind Axes 


3 " 

B 


directed along the relative wind 

normal to the wind and in a plane normal to the plane of symmetry 
normal to the wind and in the plane of symmetry, directed "down" 
origin in space 


Stability Axes 

Y s in the plane of symmetry and coincident with the projection of the 

wind in that plane 

Cj s normal to the plane of symmetry, directed along the right wing 

^ in the plane of symmetry, directed "down" 

B origin in space 


Body Axes 

X in the plane of symmetry, directed toward the nose of the aircraft 

<J normal to the plane of symmetry, directed along the right wing 

^ in the plane of symmetry, directed "down" 

B origin in space 


Earth-Surface Axes 


X c $ perpendicular to the direction of the local gravitational force, 
directed north 

u perpendicular to the direction of the local gravitational force, 

es directed east 

_ coincident with the direction of the gravitational force, directed 

toward the earth 

8 origin in space 


vii 



NAV TRADEVCEN 801A 


The following axes are right-hand orthogonal sets, with origin at the center 
of the earth. 


Earth Axes 


X E in the equatorial plane, normally passing through the zero 

meridian (Greenwich) 

)£ in the equatorial plane 

E f normal to the equatorial plane, passing through the North Pole 

£ origin in space 

Earth-Centered Inertial Axes 

X 0 in the equatorial plane, normally in the direction of vernal 

equinox 

Y g in the equatorial plane 

E 0 normal to the equatorial plane, passing through the North Pole 

£ origin in space 

The following axis system is assumed to be fixed with respect to the "fixed 
stars", and for practical purposes, its origin may be considered to be the center 
of the sun. 

Inertial Axes 

Xj parallel to X 0 , normal in the direction of the vernal equinox 

Y I parallel to Y a 

E r parallel to E a 

Z origin in space 

DIMENSIONAL UNITS 


distance - feet 
time - seconds 

angle - radians (unless otherwise noted) 
force - lbs 
moment - ft-lbs 
mass - slugs 

The angles ( X , P , Ct , /8 ,8 ) and angular rates (■/> , q. , r ) are pre¬ 
sented in units of deg and deg/sec in the time histories and some figures. However, 
in all other cases (equations, derivatives, etc. ) the standard dimensions given in 
the above table are used. Any deviation always has the proper units denoted. 
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NOTATION 

[ ] 

C) 

(A); • ( 

t 

SYMBOLS 

7* 

w 

F i F~ x $ Fy, Fj 
M t M x , M y , Mj 
>yt 

-?i Jx, 

V ; a, sir, or 

7 °* 9> r 

Ven 

*>e i f>* » f t , r* 
JL 

F 


^8 

W; W.,,W x ,-H r 


y“i s“x »s“y 
Va. i 

ft *9+’**■ 


5, m 9 

F r ,Mr 
L, 0 
*,0 


matrix 

vector (superposed bar) 

time rate of change with respect to the " i " axes 
time 


linear momentum of the body 

angular momentum of the body about the c.g. 

sum of forces on the body; body axes components 

sum of moments on the body about the c. g. ; body axes components 

mass of vehicle 

moment of inertia operator relative to c. g. ; body axes components 
inertial velocity of body axes (not including earth's orbital velocity); 
body axes components 

inertial angular ocity; body axes components 

relative veloctiy of body with respect to earth, earth-surface axes 
components 

inertial velocity of body due to earth rotation 

angular velocity of body with respect to earth; body axes components 
inertial angular velocity of the earth 
position vector 

position of earth with respect to inertial axes 
position of body with respect to inertial axes 
position of body with respect to earth axes 

velocity of wind with respect to earth, earth-surface axes 
components 

gust velocity, body axes components 

velocity of body axes with respect to air, body axes components 

angular velocity of body axes with respect to the local air, body 
axes components 

aerodynamic forces and moments 
gravitational forces and moments 

forces and moments due to thrust or reaction controls 
lift, drag 
longitude, latitude 
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2 

A 

R 0 

^ 6 , 0 
%,e s ,4 
a. 

fi 

S e , S At s r 
5, b, c 



T 

P 

M 



Q 

R 

c 


tangential velocity over the earth's surface 
heading angle of trajectory over the earth's surface 
climb or flight path angle 
altitude 

radius of earth 

normal Euler angles of body axes with respect to earth-surface axes 
special Euler angles 
angle of attack 
sideslip (angle) 

elevator, aileron, and rudder deflections 


wing area, span, and reference chord 
force coefficients, body axes 
moment coefficients, body axes 
lift and drag coefficients, wind axes 
acceleration of gravity at sea level 
atmospheric density 
speed of sound 
atmospheric temperature 
atmospheric pressure 
Mach number, V A /a, 

Reynolds number 


incompressible dynamic pressure 
pressure coefficient, — ^* rAT -- 
stagnation point heat transfer rate 


nose radius 


natural frequency (undamped) 
damping ratio 
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Body Axes 


LIST OF MATRICES 



external forces 



external moments about c. g. 



components of inertial angular velocity 



components of inertial angular velocity 



moments and products of inertial about the 
c.g. - X - “J. symmetry assumed 



components of angular velocity with respect 
to earth axes 
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W- 






components of velocity oi c. g. with respect 
to air 



components Of angular velocity with respect 
to ait 1 


M 


A 

A 


l/W 


components of velocity of the air with respect 
to the earth axes - the increment over and 
above general air mass motions (gusts, etc. ) 
- defined with respect to aircraft 



aerodynamic force coefficients 



aerodynamic moment coefficients 


Earth Axes 


M 


o 

o 

n 


angular velocity of earth with respect to 
inertial space 
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ft 
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w- 


components of the position of the c. g. with 
respect to earth axes 


Earth-Surface Axes 



components of the velocity of the c. g. with 
respect to earth axes 





the inertial velocity of the vehicle due to 
earth rotation: = Jl r* aB c*t i> 


[w> 



components of velocity of the air - general 
northward, eastward (westerly), and upward 
air mass motions (winds, etc, ) - defined 
with respect to the earth 



components of the position vector of the body 
axes with respect to earth axes 


I 
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Orientation Angles 



position coordinates of the earth-surface 
axes with respect to earth axes: A , altitude; 
X . longitude; , latitude 



successive rotations defining the orientation 
of the body axes with respect to earth-surface 
axes: ft, azimuth ("yaw"); Q , elevation 
("pitch"); 0 , bank angle ("roll") 



alternate rotations defining the orientation 
of body axes with respect to earth-surface 
axes - used for dives, spins, high angle of 
attack re-entries, etc. : ft t , wing azimuth 
or spin angle ("roll"); , nose tilt ("pitch"); 

0 5 , wing tilt ("yaw") 


Transformations 



-M MW 

transformation from earth-surface to body axes 

M 

-HWWW 

alternate transformation from earth-surface 
to body axes - used for • —* ±90* 

[**-*•] 

-M -M 

l = 1, 2, 3; J = 1. 2, 3 

M 

-MM M 

transformation from earth axes to earth- 
surface axes 

M 


transformation from earth to body axes 


-HM 

transformation from wind axes to body axes 
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M- 

[®] 

M = 



i 

0 

0 ' 

H- 

0 

0 

/ 

= 

0 

cos 0 

sin 0 

0 

/ 

O 


0 

-sin 0 

cos 0 


-7 

0 

0 


Cos 9 

0 

-sin 9 

M- 

60S ^ 

0 

sin 1> 

s 

0 

/ 

0 

0 

/ 

O 


so 1 0 

O 

cose 


-5/W 

o 

cosl> 


~C0S yr 

sin p 

o ~ 


* 60S A. 

sin X 

o 

= 

- sin yr 

cos Hr 

0 

W- 

-S/>» A 

cot X 

O 


0 

O 

/ 


0 

O 

/ 


[**] is a 90“ negative (counterclockwise) about the y -axis. 

M ■ M • Ps] have the same form as M • W , and Or] ■ 

respectively. 


M 


C3*/0 -Hn/S 0 

s infl cos /S o 

0 O / 


w- 


60S df 

o 

S/* 4 


0 

/ 

0 


-*/'* & 

o 

cos 4 
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Operators 

The general matrix 


w- 


* 

*/ 


is used to illustrate matrix operational notation. 
Cross Product Operator: 

Transpose: 

[»>]' * [*'*] 


° -f y 

£ © -X 

-y % o 


Derivative: 


w- 


y 
i. 


Inverse: 

M * V*-*\ [*•-«■] 1 * M is the unity matrix, 
rwr- hj ^ 1 for orthogonal transformations. 


o o 

/ © 

© / 

is the same function of 0 % , as [ C ] is of © 
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Combinations and Special Matrices 


M 


CM d 


«/«g 

.C#l« 
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KJ ■ [ 


O 

J1 t/tt i> j 




o' 

V 


X ctf *>' 

[^«4, 

J* 

- i> 

_ 0 

♦ MM o 
L*J 

• 

-1> 

-X */**> 


EMJ • [a*][»•«] * ^ 




-Y„ 

**• 

o 



0 tf*1> 

0 ’ 

V 

* [(’•ft).* '] [*] - -n- 

- tin 0 

O c»t1> 

0 

;Vr 
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SECTION I 
INTRODUCTION 


1. 1 BACKGROUND 

The rocket engine, with its tremendous thrust capabilities, has pro¬ 
duced a quantum jump in the possible aircraft flight speeds. Hypersonic flight 
for manned aircraft, at speeds near orbital velocity, can now be considered 
practical, and manned space vehicles are predicted for the near future. Hy¬ 
personic aircraft, with range capabilities from a few thousand miles to a few- 
rotations of the earth, are the first breed of vehicles being considered. Ex¬ 
amples are the X-15 airplane and the vehicles being developed under the Mer¬ 
cury and Dyna-Soar programs. 

The role of man in hypervelocity flight vehicles has been the object of 
considerable interest and conjecture. Traditionally, man has been an element 
in the primary airplane control loop, but recent trends place more reliance on 
automatic controls and relatively less on man for the primary control function. 
This trend has been necessary because of man's limitations and the more strin¬ 
gent control requirements, rather than because of a desire to relieve or re¬ 
place man's functions. Man's performance as an element in a control system 
is in certain ways very poor: his motor actions are power limited and slow; he 
requires an environment carefully controlled in temperature, pressure, and 
physical and chemical composition; and his tolerance to accelerations or sit¬ 
uations producing mental stresses is limited. But in other ways his perform¬ 
ance is unexcelled: he is extremely adaptable and is able to accept a multitude 
of different kinds of inputs, organize them, and interpret their meaning in terms 
of corrective action; man represents an extremely capable and efficient logical 
computer. Thus, if man's capabilities are to be efficiently employed, we must 
provide him with environmental conditions compatible with good performance 
and present him with the proper cues so that he can perform assigned tasks. 
These assigned tasks must be compatible with his limited motor capabilities. 

For man to perform effectively, he must not be overloaded with stabilization 
and control tasks which could be better performed by automatic devices; rather, 
enough of his potential must be reserved for tasks which utilize his unique de¬ 
cision making and adaptive capabilities. 

Generalities of the foregoing type are rather easily come by, but it is 
not such a simple matter to specify the requirements for a particular vehicle 
or particular situation or set of circumstances. If one is dealing with hard¬ 
ware, one can usually find a mathematical model with which to represent the 
system and which describes system performance with sufficient accuracy to pro¬ 
vide the system designer with all the information he needs to come up with a 
satisfactory design. However, man is so complex an element that such an ap¬ 
proach cannot be used. Instead, one tries to accumulate data on particular de¬ 
tails of the problem. For example, dimensional data on various pilots is used 
to design the physical layout of a cockpit and enough flexibility is built in to pro¬ 
vide for the deviation from the normal values. Another approach has been to 
try to characterize man's actions under very special circumstances by an equa¬ 
tion. Thus, work has been done on defining transfer functions for the human pilot 
in order to determine how he operates as a stabilizing device in an airplane, and 
what the characteristics of the airplane should be without the man in the loop. 
However, these are very specialized studies related to very particular circum¬ 
stances. More generally, a large body of information has been collected by 
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engineers, psychologists, and bio-medical people. Based on past experi¬ 
ence, one uses this data to try to design a compatible man-machine system. 
Then one builds such a system and modifies it to correct the errors. 

Simulators offer a useful tool in system design. Rather elaborate 
simulation facilities can often be justified because of the expense involved in 
building, operating, and testing prototype flight systems. There are three 
types of simulators. One is the research simulator which is used to obtain 
general characteristics and data not related to specific systems. The vari¬ 
able stability airplane is a good example of such a research simulator. An¬ 
other type of simulator is the design or development simulator used in the de¬ 
sign of a specific vehicle. Finally we have the training simulator which is 
used to actually train pilots or improve their proficiency in flying specific 
airplanes. All three types of simulators are needed for hypervelocity flight 
vehicles. Much work has been done in this field, particularly in connection 
with the X-15, Dyna-Soar and Mercury projects; more is needed. 

The Naval Training Device Center has sponsored, as a joint Navy and 
Air Force project, the development of a high-speed general purpose digital 
computer for application to real time simulation. The concepts were originated 
at the Moore School of Electrical Engineering, University of Pennsylvania; 
the Sylvania Electronic Systems Division of Sylvania Electric Products, Inc. 
developed the concepts into UDOFTT (Universal Digital Operational Flight 
Trainer Tool). UDOFTT consists primarily of the stored program digital com¬ 
puter, two aircraft cockpit mock-ups, and control consoles. This facility, 
located at Garden City, Long Island, is being used as a research tool for in¬ 
vestigation of various real time simulation problems. The Naval Training 
Device Center, as part of a program to exploit the capabilities of UDOFTT, 
has sponsored at the Flight Research Department of the Cornell Aeronautical 
Laboratory, Inc. a study of the equations of motion required to simulate the re¬ 
entry and landing of manned hypersonic flight vehicles of the Dyna-Soar type. 

1.2 PURPOSE 

The end objective of performing this study, as alluded to previously, is 
to determine the feasibility of using UDOFTT as a piloted simulator for hyper¬ 
sonic flight vehicles. However, the specific results of this study, as presented 
in this report, are not intended to answer this question directly, but rather are 
intended to present just the necessary facts relative to the equations of motion. 

On the surface, it might appear that six-degree-of-freedom equations 
of motion for a re-entry vehicle would not be substantially different from those 
for a supersonic fighter. However, a more careful examination of the problem 
shows this not to be the case. The new geophysical aspects entering the prob¬ 
lem, such as spheroidal earth, earth rotation, wind gradients, atmospheric 
properties at high altitudes, and gravitational gradients, all complicate the 
equations of motion. The large variations in vehicle attitude often associated 
with re-entry (e.g. , angles of attack to 90 s ) and the non-linearities associated 
with hypersonic aerodynamics complicate the aerodynamic representation. 
Furthermore, the extreme speeds and altitudes characteristic of re-entry pose 
severe requirements on computational accuracy for solving the equations. 

Recognizing that most of the complexity is introduced by the new 
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dimensions associated with near orbital and re-entry flight, it became the pur¬ 
pose of this study to investigate primarily these regions of flight. Also, in 
recognition of the complexity of the problem, it became the purpose of this 
study to find what simplifications could be made in the equations and still re¬ 
tain adequate simulation. Furthermore, since the equations were to be devel¬ 
oped for a piloted maneuvering vehicle in which the primary emphasis was on 
the piloting problem (in other words, the solution was required to be accurate 
only so far as it affected the behavior of the pilot), it was recognized that anal¬ 
ysis of the pilot's task would have an important influence on the required ac¬ 
curacy of the equations and the acceptability of the various simplifications. 

Thus, the purpose of this study is to provide a set of equations of mo¬ 
tion suitable for simulating manned hypersonic flight on. a digital computer 
such as UDOFTT. Since the simulation is concerned primarily with the piloting 
task, the equations are intended to represent the re-entry and landing portion 
of the flight regime. Because the new aspects of the problem are associated 
with the near orbital portions of this flight regime, most of the study has been 
concerned with such problems. Little work has been done on low speed prob¬ 
lems as these are already well understood. 

Also, since the configurations and missions of hypersonic vehicles are 
not at present well defined, and since there are still large areas of hypersonic 
flight about which little is known, it is not reasonable to expect that any one par¬ 
ticular set of equations will satisfy all requirements. For this reason, it is 
important to summarize the basis on which the recommended equations have been 
formulated and to summarize methods by which they can be extended, simplified, 
or altered to account for future changes in the circumstances governing their 
use. To this intent, sections on geophysics and aerodynamics, a complete der¬ 
ivation of the equations of motion, and appendices on extension of the equations 
have been included in this report. 

1.3 APPROACH 

The primary tasks in this study have been to derive equations of motion 
for a maneuvering re-entry vehicle, program these equations for solutions on 
an IBM 704 digital computer, and compute solutions of these equations to de¬ 
termine the effects of various parameters in the equations and the effects of 
making various simplifying assumptions in the equations. Early in the study, 
visits were made to several aircraft companies and NASA to discuss pertinent 
work being done on re-entry, with respect to both the simulation problem and 
the characteristics of manned hypersonic aircraft. The information so obtained 
was of real value in guiding the direction of study. The one outstanding finding 
was that (though several groups had attempted or were starting programs to 
solve the complete equations of motion for hypersonic flight) no such program 
existed in a developed state*. Numerous re-entry programs existed, but all 


* Since completion of this study, the recently completed work of Brown, Brulle 
and Griffin was uncovered and their basic development of the equations of motion 
parallels closely the present work. While Brown et al. develop a generalized 
computer program in great detail, the present work was directed toward sim- 
lifying the equations to facilitate their solution on a real time digital computer. 
Brown, R.C., R.V. Brulle and G. D. Griffin: "Six-Degree-of-Freedom Flight- 
Path Study Generalized Computer Program - Part I, Problem Formulation: 
WADD Technical Report 60-781, November I960. 
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involved substantial simplifying assumptions, primarily in the degrees of 
freedom considered or in the representation of the earth's shape, rotation, and 
gravitational force. It was apparent from our own initial efforts and from dis¬ 
cussion with other people in the field that the task of developing, programming, 
and checking out complete equations of motion for a maneuvering hypersonic 
aircraft was no mean feat and that finding acceptable simplifications in the 
equations was of paramount importance. 

Accordingly, equations of motion were derived which represented the 
conditions of near orbital and re-entry flight as completely as possible, but 
with certain notable exceptions. The earth was assumed spherical and.no thrust 
or mass variations were included in the equations. Earth oblateness was ex¬ 
cluded (both the geometric effect and the gravitational field effect) because it 
complicated the process of checking out the computer program and analyzing 
the results, and because sufficient work has been done by other investigators 
to show when oblateness has to be included. Variable mass and thrust were ex¬ 
cluded arbitrarily because the vehicle was assumed to be a boost glider that was 
unpowered during re-entry, though this would not exclude the concept of a landing 
engine. 


A typical re-entry vehicle configuration was selected for study - one with 
a low aspect ratio modified delta wing with vertical tails mounted at the wing tips 
and with fuselage above the wing. Estimates of lift, drag, and inertial charac¬ 
teristics were obtained and an idealized control system was selected to provide 
control of angle of attack, sideslip, and bank angle. The rapid dynamic response 
so obtained would be similar to the response that could be realized if the vehicle 
had an adaptive control system. Thus, though no aerodynamic moments as such 
were included, a useful and reasonable simulation of the over-all aerodynamic 
characteristics was obtained. 

The principal area of interest for this study was defined as the re-entry 
and landing of hypersonic flight vehicles; more precisely, flight below 300,000 
feet of altitude and less than 26, 000 ft/sec (approximate circular orbital ve¬ 
locity). Accordingly, the geophysical models of interest are those for the 
earth's shape, the gravitational field, and the atmospheric density, tempera¬ 
ture, winds, and turbulence (gusts). There is adequate data to formulate ten¬ 
tative models, within the defined limits, for all but atmospheric winds and 
turbulence. Existent data on winds at high altitudes is mostly for the northern 
hemisphere, particularly over the United States; almost no data is available 
on turbulence at high altitudes. Most of the atmospheric data is based on mean 
values and the data on temporal variations, particularly short time variations, 
is fragmentary. In Section II of this report, the present status of our knowledge 
of the pertinent geophysical characteristics is summarized, and mean or typ¬ 
ical.models are formulated for these characteristics. The air has been treated 
in this study as a normal continuum since free molecular flow effects do not 
take on significance until altitudes above 300, 000 feet are reached (References 
1.4, 1.5, 2.5). 

A complete formulation of the aerodynamic characteristics for hyper¬ 
sonic vehicles cannot be made at the present time; adequate theory and data 
are not available. Even an attempt to discuss the available material is far 
beyond the scope of this investigation. However, enough effort was devoted to 
aerodynamics so as to define the over-all problem. Re-entry aerodynamics 
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is complex, more so than that for low speeds because of the new facets in¬ 
volved. There are three basic aspects: (1) lift and drag, the performance 
parameters which govern the vehicle's trajectory and provide path control, 
and the vehicle's deceleration which must be maintained within limits toler¬ 
able to the pilot; (2) aerodynamic heating, essentially a local flow problem 
which poses limits on the allowable flight regime; (3) stability and control 
characteristics, involving complex over-all flow problems, which can be ex¬ 
pected to lead to stringent requirements for automatic control and limits on 
the allowable flight regime. The aerodynamic problem is compounded by the 
following aspect. We are dealing with a many faceted air: subsonic, tran¬ 
sonic, supersonic, and hypersonic flow; ionization and dissociation, both lam¬ 
inar and turbulent boundary layers, and free molecular flow at the higher al¬ 
titudes. We are dealing with a strongly non-linear problem: the basic non¬ 
linear nature of hypersonic flow and the high angles of attack contemplated in 
re-entry will cause severe non-linear stability and control problems. Fur¬ 
thermore, in hypersonic flight the problems associated with performance, 
heating, and stability and control are intimately related and cannot be treated 
independently, and their solutions are not mutually compatible. If the past 
can be used as an indicator, vehicle configurations will be established prin¬ 
cipally by performance and heating considerations, leaving stability and con¬ 
trol to make the most of a poor situation. On this basis, we can expect strin¬ 
gent requirements for automatic control; should the automatic control system 
fail, it will be enough to just return man safely to earth - if even that is possible. 

Much of the available data on projected manned hypersonic vehicles is 
classified, but the general aerodynamic characteristics can be formulated based 
on the available unclassified theory and experimental data. The vehicle's basic 
over-all characteristics at hypersonic speeds (lift, drag, stability derivatives) 
are similar to those at lower speeds, but non-linear effects are abnormally 
strong. We have attempted in this study to formulate a guide for suitably rep¬ 
resenting the aerodynamic data in the equations of motion, based on the re¬ 
quirements of the simulation task. 

As for the control problem, it will be severe since the uncontrolled ve¬ 
hicle can be expected to have poor stability and control characteristics and 
both aerodynamic and reaction controls will be necessary. An adaptive control 
approach will very likely be used because of the extremes involved in the con¬ 
trol problem. Achieving adequate reliability will be difficult. Performance and 
aerodynamic heating will enter the control problem in the form of energy man¬ 
agement and temperature control systems and pilot's displays. A unified ap¬ 
proach will be needed for the whole control problem, including flight controls, 
energy management system (including propulsion), temperature control system, 
navigation and guidance system, and the pilot with his displays and controls. In 
this study, we have not considered all these problems in detail. Rather, we 
have tried to formulate the equations of motion so that the necessary data (var¬ 
iables and quantities) would be available for the pilot's displays and the various 
systems. 

The derived equations were programmed on an IBM 704 digital computer, 
and a fourth order Runge-Kutta numerical integration scheme was used to solve 
the equation. Work was done on a more sophisticated integration technique util¬ 
izing a "predictor-corrector" integration formula which automatically varied 
the integration time interval. However, there was insufficient time to develop 
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the numerical error criteria required to adjust the time interval, and the 
straight Runge-Kutta integration technique was used for all the solutions ob¬ 
tained in this study. Considerable difficulty was encountered in checking the 
equations as programmed on the IBM 704. The normal difficulties (elimi¬ 
nating mistakes in the equations and the program) were compounded, by the 
difficulty in obtaining a satisfactory numerical integration technique. It was 
estimated that if the Runge-Kutta method was used and the control system par¬ 
ameters were selected to give good handling qualities, it would take 10 hours 
of computing time to solve for one re-entry - a prohibitively long time. This 
long computing time evolved because the integration time interval had to be 
short enough to compute the fastest mode of motion present in the equations 
(the longitudinal and lateral short period motions which for good handling qual¬ 
ities would have 0.5 cps natural frequencies). If such solutions had been com¬ 
puted, it was not certain that they would have been correct because of the 
round-off errors in the integration of the variables describing the trajectory 
of the vehicle. The scope of the work did not permit an extensive study to 
solve these difficulties so a useful expedient was adopted. It was reasoned 
that the combined dynamics of airplane and pilot, as the pilot tried to control 
his trajectory in re-entry would have a much lower effective natural frequency 
than that of the airplane alone. It was determined that a . 05 cps natural fre¬ 
quency was typical of the control-fixed longitudinal short period mode for a re¬ 
entry vehicle at 250,000 feet of altitude (Reference 1.6). Accordingly, the con¬ 
trol system parameters were selected to provide a 0.7 damping ratio (g ) and 
.05 cps natural frequency ( co n ) in the angle of attack, sideslip, and bank angle 
modes. By this reduction in short period natural frequency from 0. 5 cps to 
.05 cps, the computing time for a re-entry was reduced to about one hour. A 
further reduction to about 20 minutes was effected by the following computational 
expedient. For those portions of the re-entry where no short period response 
were being excited by command inputs and the control system was simply reg¬ 
ulating, an integration time interval was selected just small enough to provide 
a stable computation. Then when command inputs were inserted, simulating 
positive pilot control action, the time interval was decreased to obtain an ac¬ 
curate computation of the short period motions. In this manner it was possible, 
using the Runge-Kutta integration technique, to obtain a reasonably short com¬ 
puting time for re-entries and at the same time compute accurately both short 
period and long period or trajectory motions. 

Considerable thought was given as to what solutions should be computed. 
It was clear that because of the purpose of this study, it was not necessary to 
hold the re-entry to any precise limits. Any solution which approximated re¬ 
entry conditions would satisfy, but on the other hand it was necessary that the 
solutions encompass all the elements that would actually be encountered (e.g. , 
non-equatorial flight, skips during re-entry, both lateral and longitudinal 
maneuvers, excitation of both short and long period dynamic modes and their 
interaction). It was also necessary to obtain some measure of significance of 
the various effects in relation to the piloting problem, particularly to assess 
the quantitative significance of simplifications in the equations of motion. Ac¬ 
cordingly, three basic solutions or cases were se*ected for study: (1) an equa¬ 
torial re-entry at constant angle of attack, (2) a similar re-entry but aligned 
down the Atlantic Missile Range, and (3) a re-entry initially like (2), but in 
which a 45* banked turn was initiated at the bottom of the first skip and held 
for the duration of the solution. Perturbations about Case (1). were made in 
which the angle of attack (oc), the bank angle ( 0 ), and sideslip ( j3 ) were 
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varied. These perturbation solutions provided a measure of the effective¬ 
ness and sensitivity of the pilot's control over the vehicle's trajectory. The 
effect of short period dynamics was determined by effecting the banked turn 
and the perturbations in at , /3 , and 0 by step command inputs, all at the 
bottom of the first skip. 

It would have been desirable to compute additional types of solutions 
(particularly re-entries from the recovery ceiling, re-entries for different 
Ofc's and L/D's, and re-entries with winds and gusts), but shortage of time 
precluded such solutions. 

The programmed equations included certain initial simplifications: no 
earth oblateness, no thrust, constant mass. The effects of oblateness are de¬ 
terminable from the literature, and their pertinence to the simulation problem 
is discussed in this report (Section 6.2). Thrust and variable mass are out¬ 
side the scope of this study, but methods for their inclusion in the equations of 
motion are discussed in Appendix B. Study of simplifications in the equations 
of motion, other than these two initial ones, has been effected by computing 
solutions on the IBM 704 with the simplified equations and comparing these so¬ 
lutions to those obtained with the complete equations. It was not possible to 
study the effect of each individual term in the equations - far too lengthy a 
process. Instead, the effect of various groups of terms with some physical 
significance were examined (e. g. , the effect of earth rotation in the equations 
fordl and £ ); and also, total effects of certain quantities were examined 
(e.g. , the effect of neglecting earth rotation entirely). It was necessary to 
assess the effect of a simplification from two standpoints: first, from the ef¬ 
fect on the over-all trajectory; and second, from the effect on the short per¬ 
iod motions. Also, it was necessary to consider interaction of the two effects. 
Evaluation of whether a simplification was allowable or not was primarily based 
on how it affected the piloting task. Here, the perturbation solutions (Case (3) 
with variations in QL , 0 , and ) helped to give some quantitative measure of 
the importance of a change in the solution; but at the same time, experience 
and knowledge of airplane handling qualities and the piloting task were also 
important. 

Recognizing that many problems have not received adequate evaluation 
in this study and that concepts change rapidly, we do not consider the equations 
developed in this report as necessarily optimum nor have we attempted to be 
completely explicit in the formulation. Accordingly, we have attempted to de¬ 
note the shortcomings of the equations, to point out possible methods for im¬ 
provement, and to show in what areas there is insufficient knowledge for defin¬ 
itive formulation. 
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SECTION II 
GEOPHYSICS 


2.1 INTRODUCTION 

Since any geophysical property of the earth may have many different 
mathematical models of varying complexity depending upon application, the 
properties which affect atmospheric re-entry flight have been studied with the 
objectives of this project as a guide; that is, the formulation of re-entry ve¬ 
hicle equations of motion suitable for use with a fixed-base flight simulator 
driven by a digital computer. It is for the purpose of defining adequate models 
of the pertinent geophysical properties that this section is included. 

All geophysical elements usually modeled for the study of normal at¬ 
mospheric flight must be considered; but, in addition, some of these concepts 
must be extended and other completely new concepts added. In other words, 
the operation of a re-entry vehicle within the normal flight range is just one 
segment of a wider range of operation which begins with orbital conditions and 
ultimately ends with low subsonic velocities near the surface of the earth. In 
the paragraphs which follow, a qualitative discussion of various phases of re¬ 
entry is given. 

In a near-orbital condition (defined here as an orbit whose lifetime is 
only a few revolutions) the major force acting on a re-entry vehicle is the 
gravitational force which is essentially balanced by the so-called centrifugal 
force due to the vehicle's velocity. Present also are the very small aerodynam¬ 
ic forces drag and lift. It is the drag force which causes inevitable orbital de¬ 
cay while the lift force produces "long" and "short" period oscillations analo¬ 
gous to the usual low altitude longitudinal flight oscillations. The analysis of 
Reference 1.6 shows that the period of the "long" period oscillation approaches 
orbital period and, in variance with its name, the period of the "short" period 
oscillation approaches infinity as altitude increases. As small as the aerody¬ 
namic forces are, however, they completely mask any effects of lunar and solar 
gravitational fields, the earth's magnetic field, radiation pressures, atomic 
particle and micrometeorite collisions, and other lesser effects. 

The gravity and drag forces are probably more important in this near¬ 
orbital condition. Accurate description of the vehicle's orbital lifetime, a 
function of drag, is dependent upon adequate knowledge of air density and ac¬ 
curate determination of vehicle velocity with respect to the air. The former 
requires a good model of density while the latter requires that the atmosphere 
rotate with the earth. This means that the earth's angular velocity should be 
included with the equations of motion. There are smaller periodic and secular 
(non-periodic) drag effects, but they are not discussed here. 

The largest single external force acting on an orbiting vehicle is gravity. 
As such, it affects all flight near the earth. With regard to an orbital vehicle, 
the earth's true gravitational field produces some small perturbation effects 
due to the oblateness in addition to the usual radially directed force. The per¬ 
turbations are periodic and secular with secular motions such as regression of 
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nodes* and precession of the line of apsides** being most important. Further 
information can be obtained from Reference 2.1. As a consequence of their 
magnitudes, the gravitationally-induced periodic effects are negligible and the 
secular effects become important only for an extended stay in orbit. Since an 
extended stay in orbit is not at all likely for;a re-entry vehicle, the gravita¬ 
tional effects of oblateness are probably not important in this study. 

The most important effect of the earth's oblateness is due to the geometry 
involved. The path of an orbiting vehicle near the earth or one re-entering the at¬ 
mosphere is affected because the oblateness causes periodic changes in atmos¬ 
pheric density for the orbiting vehicle and similar but not necessarily periodic 
changes in density for a re-entering vehicle. Earth oblateness also produces, 
in analysis, a vehicle location problem because pseudo-spherical*** coordinates 
(7i , \ ) such as derived in Section III no longer give the vehicle's location 

with respect to the earth. Oblateness can cause latitude errors of nearly 12 
miles on the earth's surface and differences of about 6. 5 miles in altitude when 
compared with a mean spherical earth. 

As the vehicle penetrates deeper into the sensible atmosphere, the phase 
actually called re-entry begins. During the re-entry the most severe aerody¬ 
namic effects, deceleration and heating, occur. Drag is indirectly responsible 
for both, in the sense that deceleration is proportional to drag and the vehicle's 
kinetic energy is, through the mechanism of drag, converted to heat. Lift force 
is an important moderating factor during re-entry because it can be controlled to 
decrease decelerations and heating rates. However, lift can, and usually does, 
cause some "skipping" which is the "long" period oscillation mentioned previous¬ 
ly, and an increase of total heat input. In addition, lift provides the capability 
for correcting certain miscalculations or errors made previously during re-entry. 
It becomes apparent that geophysical properties must be adequately modeled for 
the re-entry phase, especially the atmospheric density and temperature models 
which are very necessary for computation of aerodynamic forces. Earth rota¬ 
tion, too, with the corresponding atmospheric rotation, grows increasingly more 
important as the vehicle's inertial velocity decreases in re-entry. 

The final phase of a re-entry might be described as "normal flight" be¬ 
cause the vehicle is operating at velocities and altitudes common to ordinary high 
performance aircraft. Here, restrictions are relaxed on requirements of some 
of the geophysical properties; for example, the gravitational acceleration could 
be taken as constant and the earth's shape can be considered spherical, for at 
this velocity and altitude the vehicle's range is short. However, it may be de¬ 
sirable to evaluate the effect of winds and gusts upon the vehicle. The large scale 
winds probably would cause deviations from the desired path while gusts could 


* The line of nodes is the line which connects points of intersection of the orbit 
and the equatorial plane. 

** The line of apsides is the major axis of the elliptical orbit, the line connecting 
the apogee and perigee which are respectively the farthest and the nearest points 
from the center of the gravitational field. 

*** The variables ( , p , X ) are called pseudo-spherical because the usual 

spherical coordinates are equivalent to ( , 90®- p , ) in the derivation of 

Section III. 
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cause difficult and dangerous handling problems for the pilot of a re-entry ve¬ 
hicle, especially if the vehicle is marginally stable. 

Any set of geophysical properties chosen for equations of motion of a 
re-entry vehicle would be consistent in themselves and therefore '•correct". 

For a flight simulator, however, the geophysical modeling should be made as 
nearly "real" as possible. The simulator models should provide, simply and 
adequately, the realism needed throughout the various phases of re-entry flight. 
The next sections provide a description of the various geophysical properties of 
the earth. The problem will be divided into two major categories: first, the 
motions, size and shape, and gravitational field of the earth; and second, the 
earth's atmosphere. 

2. 2 THE MOTIONS, SHAPE, AND GRAVITATIONAL FIELD OF THE EARTH 

In order that accurate equations of motion for a vehicle may be written, 
the equations must be referred to an inertial frame. For orbital or re-entry 
studies a convenient reference frame would be centered in the earth and possi¬ 
bly rotating with the earth. However, this earth-centered frame is not an 
inertial frame because the earth is not fixed in space. It is possible to refer 
all motions of a vehicle to such an earth-centered frame provided the motions of 
the earth are known. 

The two principal motions of the earth which might significantly affect 
the equations of motion are: 

1. rotation about its axis, and 

2. revolution about the sun. 

These motions will be discussed in detail later with emphasiB placed upon appli¬ 
cability of the motion to this study. Additional information is given by Parvin, 
Reference 2.2. 

It is necessary to know the shape or figure of the earth so that the posi¬ 
tion of an orbiting or re-entry vehicle may be accurately described. Generally, 
this is a more difficult situation to resolve mathematically than the motions of 
the earth. Altitude can be obtained in terms of the vehicle's spherical coordi¬ 
nates for an earth-centered frame by knowing the magnitude of the radius vector 
and the radius of the earth at some point on its surface. It is the model of the 
earth's size and shape which defines the radius of the earth. Altitude is used in 
preference to the radius vector magnitude because the earth's surface provides 
a convenient reference for altitude determination, the concept of altitude has 
been in general use for aeronautical work since the earliest days of flight, and 
atmospheric properties, as described in Section 2. 3 for example, are presented 
in terms of altitude. Unfortunately, altitude must be described mathematically 
in a different way for each model of the shape of the earth, although the defini¬ 
tion of altitude remains the same in each case. Altitude is generally defined as 
the height of a vehicle above some reference spheroid, measured normal to the 
local surface of the spheroid. 

Latitude and longitude of the vehicle are mathematically described, like 
altitude, differently for each earth-shape model. In general, there exists more 
than one definition of latitude, but usually only one definition of longitude. The 
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latter may be accounted for by the assumed axially symmetric figure for the 
earth. Longitude is defined as the angle between the projection of the radius 
vector on the. equatorial plane and some reference direction in that plane which 
is usally the direction of the line from the earth's center through Greenwich, 
England. Two definitions of latitude are given here. Geocentric latitude, as 
its name may imply, is the angular measure from the equatorial plane to the 
radius vector of some point using the earth center as vertex for the angle. 
Geodetic latitude is the angle between the equatorial plane and the line to the 
given point which is normal to the reference spheroid. A more complete an¬ 
alysis of the vehicle position problem is given in Appendix C. 

In many analyses of flight the accelerat jt\ of gravity can be taken as a 
constant without serious restriction; however, to study the motions of a vehicle 
during atmospheric re-entry and during its orbital and boost phases, a more 
accurate description of this acceleration is undoubtedly necessary. The possi¬ 
bilities include a simple inverse-square variation and the more complicated but 
more nearly correct oblate spheroidal variation which has a non-radial compo¬ 
nent. An adequate description of gravitational acceleration for use with a re¬ 
entry flight simulator has been an objective of this study and is reported here. 

2.2.1 Earth Motions 


The earth's rotation is clearly the most significant of all its motions so 
far as equations of motion for a hypersonic re-entry are concerned. Rotation 
about the earth's polar axis is a full 360° each day in a counterclockwise sense 
when viewed from above the North Pole. With respect to the so-called "fixed" 
stars, a 360° rotation is described as a sidereal day which is 23 hours, 56 
minutes, 4.09 seconds in length. With respect to the sun, however, the full 
rotation of the earth is more than 360* in space. The solar day is, therefore, 
longer than the sidereal day by the time required for the extra rotation: about 
4 minutes, but not exactly so because each solar day is of different duration. 

This effect is due to the earth's orbital motion about the sun. The durations 
mentioned above are in terms of mean solar time, which is identical with ordi¬ 
nary civil time except for its reference. Mean solar time will be used in all 
discussions henceforth. 

One complete revolution of the earth in its orbit about the sun, with re*- 
spect to inertial space, requires 365.256, 360,42 mean solar days* and is called 
a sidereal year. The orbit is slightly elliptical with an eccentricity of 0.0167272, 
and the earth's revolution is counterclockwise when viewed from above the North 
Pole. It is because the earth's orbital motion and its rotation are in the same 
direction that the solar day is longer than the sidereal day. In addition, the 
ellipticity of the orbit means that the earth moves along in its path a different 
angular distance each solar day, thus accounting for the variable length of the 
solar day mentioned in the preceding paragraph. The aphelion, or the greatest 
distance from the sun, is 1.519 x 10*1 meters while the perihelion, the least 
distance, is 1.469 x 10** meters. A situation of importance is that throughout 
the earth's orbital motion the sun's gravity will attract the earth and any body 
near it (e.g. , a re-entry vehicle) with nearly equal forces. 


* The mean solar day, which has a constant length, represents the average of 
all true solar days during one revolution of the earth about the sun. 
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An elementry mathematical study of the earth's motions will show that 
only two motions, the earth's rotation around its polar axis and its revolution 
about the sun, have any appreciable effects that need be considered in the study 
of a re-entry vehicle. That is, precession and nutation of the polar axis, etc. 
can be neglected. The preceding paragraph concludes with a verbal argument 
for neglecting the translational motions of the earth, and Section 3.2.1 contains 
a mathematical description of the same argument. The result, then, is that 
the earth's sidereal rotation is the only motion that need be considered. Since 
the earth completes one rotation in space each 23 hours, 56 minutes, 4.09 sec¬ 
onds, its angular velocity is 

= 0.004178075 deg/sec. 

The earth's angular velocity is very slowly decreasing with time, but this secu¬ 
lar variation of J2 is of no concern here. Thus, by neglecting the earth's rev¬ 
olution about the sun and using the sidereal rotation, a non-rotating reference 
frame whose origin is fixed at the center of the earth may be used as an inertial 
frame. 


The earth's angular velocity causes a vehicle moving northward to appear 
to veer eastward unless it is constrained to the earth's surface. This apparent 
acceleration, that must be considered if one attempts to use the rotating earth 
as an inertial frame, is called the Coriolis acceleration. In vector notation, the 
acceleration is stated as the cross product 

A c * - 2J 2 x V £ (2.D 

in which Vg is the vector velocity with respect to the rotating earth. This 
Coriolis acceleration, a component of the total acceleration, is a direct conse¬ 
quence of twice differentiating the vehicle's position vector, , when the po¬ 

sition vector is expressed with components in an earth-centered frame rotating 
at angular velocity-A in inertial space. 

Model of Earth's Motion: The earth's angular velocity given above was 
rounded off fo^the value 


= 0.004178 deg/sec 

for obtaining the solutions of Sections V and VI. 

2.2.2 The Size and Shape of the Earth 


Neglecting surface irregularities of large and small scale, the earth is 
oblate spheroidal; that is, its polar radius (coincident with the axis of rotation) 
is smaller than its equatorial radius. It is possible to mathematically define the 
spheroid as an ellipsoid of revolution symmetric about the polar axis. In terms 
of the equatorial radius ( A ), flattening (£ ), and geocentric latitude ), the 
radius of the earth is: 





A*U-€>* 


( 2 . 2 ) 


for which adaptation of Kaula's data (Reference 2. 3) gives 
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C g a>- b - 1 

A 298.24 

& = 6, 378, 163 meters = 20, 925, 732 feet. 

For comparison with the equatorial radius, the polar radius (£ ) can be computed 
using O. and € : 

b = 6, 356, 777 meters = 20, 855, 567 feet. 

This oblate spheroidal representation of the earth's size and shape is adequate 
(or more than adequate) for nearly all applications. Equation 2.2 can be simpli¬ 
fied, if desired, by obtaining its power series expansion in terms of flattening and 
using only a few terms of the expansion. 

If a less precise model of the earth's size and shape will suffice, the earth 
can be considered spherical. It is best, then, to choose some mean radius for the 
earth in terms of the equatorial and polar radii given above. Many mean radii can 
be defined according to various criteria; however, a suitable mean is that which 
makes the volumes of sphere and spheroid equal: 

* OsO-G)* ** a*(/--$-6) (2.3) 

^ 20,902,343 feet 

Since the earth's geometric shape affects aerodynamic forces on a vehicle (indi¬ 
rectly through altitude dependent atmospheric properties such as density and temp¬ 
erature), in a particular case one might prefer a "mean" radius chosen to provide 
a better value of altitude (H ) and geocentric radius (J^ # •Xi k ^’ 4 ) than that given 
by Eq. 2. 3. This would be especially true for re-entries confined to equatorial 
or polar regions where the earth radius most differs from the mean of Eq. 2. 3. 

A more precise description of the shape of the earth has been obtained re¬ 
cently from studies of the motions of artificial satellites. Although gravity mea¬ 
surements at the earth's surface have for years indicated the presence of large- 
scale undulations of the surface, data was sufficiently scarce to prevent good es¬ 
timates of the magnitude and extent of the undulations. The artificial satellites 
have produced good gravitational data which can be analyzed to obtain a descrip¬ 
tion of the earth's shape. This data has verified that the equator is not circular, 
but slightly elliptical, making the earth more nearly a triaxial ellipsoid. In ad¬ 
dition, the "pear-shaped" figure of the earth was "found" by analyzing satellite 
tracking data. In general, these more precise models of the earth are unneces¬ 
sary for work involving equations of motion of re-entry vehicles. 

Model of Size and Shape: The solutions of the equations of motion pre¬ 
sented in this report have been computed using a spherical earth model whose 
radius is: 

= 20, 860,000 feet. 

This value is somewhat less than that indicated by Eq. 2.3, but this does not af¬ 
fect the significance of the solutions. 
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2.2.3 The Gravitational Field 


The earth's gravitational field is in reality an extremely complex thing 
which depends entirely on mass distribution. It is possible, however, to ex¬ 
press mathematically the gravitational potential to any degree of precision with 
a summation of spherical harmonics for the potential: 

U = 2 — *77 (u**) + 8* (2.4) 

■Hmo yn.%0 

in which any term, , is of degree n and order m, and (.W?) is the 

associated Legendre function. Recent studies of artificial satellite motions and 
other geodetic data have resulted in evaluation of this potential function to terms 
of m = n = 8 (see Kaula, Reference 2.3). Gravitational acceleration is then ob¬ 
tained as the gradient of the potential. 

Certainly, so complete a description of the gravitational field as is given 
above is not needed for the equations of motion developed in Section III. The po¬ 
tential function can be considerably simplified by neglecting the dependence on 
longitude (A ), a reasonable assumption for the earth. In this case, the poten¬ 
tial function, in an alternate notation, is 

V t. <*•’> 

in which G is the universal gravitational constant, is the mass of the earth, 
is the magnitude of the position vector a is the equatorial radius, and the 
J's are constants which specify the magnitude of the deviations from spherical 
(i.e., for n H 2). Kaula (Reference 2.3) states the following values: 

A = 6,378, 163 meters = 20,925,732 feet 

+3.98602 x 10 14 m 3 /sec 2 = 1.40765 x 10 16 ft 3 /sec 2 

T‘° 

J z = -1082.61 x 10‘ 6 
X = +2.05 x 10- 6 * 

T = +1.43 x 10- 6 

V 

Actually, he gives values for Jg- and but the corresponding errors are large 
relative to the coefficients themselves and, therefore, the coefficients are not 
listed here. 

Probably the most widely used description for the earth's gravitational 
field is the familiar inverse-square acceleration whose potential is given by 

JJ - (2.6) 

___ 

* This is the term which results from the earth's "pear shaped" figure. The 
small or stemmed end of the pear is at the North Pole. 
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which is simply the zero order term (n = 0) of the two more precise formula¬ 
tions given by Equations 2.4 and 2.5. The gradient of this potential gives the 
following expression with its equivalent in terms of acceleration of gravity at 
mean sea level: 



This equation clearly indicates that 6?*«« # and are equal in magnitude. Thus, 

by using the value for ^'M^given above and the value for from Eq. 2.3, the 
sea level gravitational acceleration can be evaluated: 

9 0 = = 32. 318 ft/sec 1 2 3 4 5 6 (2.8) 

% 

Of course, the expression can be used to obtain ^ if is available. 

Usually an inverse-square law for gravity is modeled in conjunction with 
a spherically shaped earth. Because there is a possibility that an inverse-square 
law would be used with an oblate spheroidal shape for the earth, it is important 
to specify what earth radius should be used to find g Since the J z term of Eq. 
2.5 is the predominant deviation from the spherical field, the other higher order 
terms can be ignored; therefore, the radius at which the second harmonic (J z ) 
term is zero ( 1) ±35. 3*) should be used with Eq. 2.8 above. 

Model of the Gravitational Acceleration: For the solutions presented in 
Sections V and VI of this report, an inverse-square gravitational acceleration 
(consistent with the spherical shape) is used. The sea level gravity is specified: 

g o = 32. 17 ft/sec 2 

a value somewhat lower than given by Eq. 2. 8. 

2. 3 THE ATMOSPHERE 

This section, which describes the atmosphere and models of some of its 
properties, begins with a discussion of general properties and atmospheric shells; 
then narrows considerably to physical properties of the atmosphere which affect 
the description of a re-entry flight, the temperature, density, winds and gusts. 

2. 3. 1 General Discussion 
--— 

The atmosphere is a very complex body of matter which can be analyzed 
(Reference 2.4) according to its various: 

1. Physical qualities (density, pressure, temperature), 

2. Chemical qualities (composition, dissociation, molecular weight), 

3. Electrical and nuclear qualities (ionization, aurora), 

4. Motions (jet stream, turbulence), 

5. Biological aspects (physiological effects upon man), and 

6. Technical aspects (aerodynamic forces and heating). 
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( 


t 


The physical qualities, motions, and the technical aspects are most important 
from the point of view of the equations of motion for a re-entry vehicle. How¬ 
ever, in the real situation, each of the six points above can be important. 
Pressure, density, and temperature are variable with altitude; and, in fact, 
temperature is commonly used as a basis for dividing the atmosphere into 
nearly concentric, approximately spherical, contiguous shells having indis¬ 
tinct boundaries (called pauses*) at various altitudes. The usual division of 
the atmosphere into shells, with a short description of some corresponding at¬ 
mospheric motions, temperature gradients, and technical aspects, follows. 

Troposphere: Winds and turbulence are prevalent; jet streams 

exist near the tropopause which is at about 45, 000 
feet. The temperature gradient is negative. This 
is the region of greatest air density. Air breathing 
engines are used and aerodynamic control of flight 
vehicles is not only possible but advantageous. 

Stratosphere: In the lower regions jet streams exist, but the 

upper regions are calm. There is no weather. 

The stratopause is placed at about 100, 000 feet, 
which also represents an approximate practical 
limit for ordinary air breathing engines. In the 
upper regions of the stratosphere, aerodynamic 
controls tend to be ineffective, which necessitates 
the use of reaction controls. The temperature 
gradient is approximately zero. 

From this point the division is not nearly so well defined. Reference 2.5 sum¬ 
marizes many different approaches to further shell separation which depend upon 
temperature, composition, ionization, etc. Here, the breakdown will follow that 
given in Reference 2.4, Table 3 - 7A, which is a breakdown based upon a com¬ 
posite of properties: 

Chemosphere: In the upper region, which has a negative temperature 

gradient, there is some turbulence and ionization 
begins. The lower region has a positive temperature 
gradient. The chemosphere is aerodynamically a 
slip flow region in which high glide-vehicle velocities 
are feasible. The region has high-velocity geostrophic 
winds, but very low density. The chemopause is 
at 260, 000 feet. 

Ionosphere: The atmosphere is increasingly ionized with altitude 

in this region due to solar ultraviolet radiation. 

There is probably some slight turbulence with in¬ 
creasing wind velocities but, like the chemosphere, 
the density is so small that the net effect on a re-entry 
vehicle is negligible. This is a region of free molec¬ 
ular flow. Temporary earth satellites are possible. 


* Pause, with the prefixed root word, is used to indicate the upper boundary of 
a region (e. g. , tropopause is the boundary between the troposphere and stratosphere). 
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The temperature gradient is positive. The iono- 
pause is at approximately 250 miles. 

These regions were essentially unexplored and rel¬ 
atively little was known about them before the era of 
earth satellites. They are what could be called re¬ 
gions of near-space flight although there is some re¬ 
maining atmosphere of ionized particles, and radia¬ 
tion belts are present. Satellite lifetimes are indef¬ 
inite. These regions are, however, of no interest 
in the usual re-entry problem. 

The above description should provide a rudimentary understanding of the environ¬ 
ment for near orbital and re-entry flights. 

2.3.2 Density and Temperature 

Since density and temperature are physically interdependent, this section 
will discuss both. The density is an atmospheric property which is of great im¬ 
portance in flight. For a given velocity, the aerodynamic forces and moments 
acting on the vehicle are proportional to the density. This explains the relative 
ineffectiveness of aerodynamic control during the first portions of a re-entry: 
even with the vehicle's great speed, the density is too small to provide the neces¬ 
sary control. Physically, the atmospheric temperature has little effect on a re¬ 
entry vehicle; but from the analysis point of view, temperature is needed so that 
speed of sound and hence Mach number can be calculated for use with the vehicle 
aerodynamic representation. 

Various methods are used to obtain the atmospheric density and temper¬ 
ature data depending upon the altitude range of interest. In the lower regions of 
the atmosphere (the troposphere and most of the stratosphere) where "air- 
breathing" aircraft can be used, the data is generally very good. Both mean 
properties and variations can be adequately measured. In the very high atmo¬ 
sphere (the ionosphere and above) earth satellites are used for determining at¬ 
mospheric properties. Although satellite data is difficult to analyze, it pro¬ 
duces good results not only for mean values of the properties but also for vari¬ 
ations of these properties. The most difficult altitude range (roughly the entire 
chemosphere) lies between these two extremes where neither ordinary aircraft 
nor satellite can operate. Here rockets, balloons and other very ingenious 
methods are used but there has been relatively little data obtained. However, 
this scant data for the chemosphere can be used to help interpolate between the 
regions on either side whose properties are more precisely known. 

The interrelationship between the measured quantities (usually pressure, 
density and speed of sound) and the desired quantities requires knowledge of the 
mean molecular weight, ratio of specific heats, and the universal gas constant. 

It is assumed that the universal gas constant applies everywhere. The ratio of 
specific heats is very nearly, constant to about 90 kilometers (« 300,000 ft); 
above this altitude, the speed of sound tends to lose significance and little is 
known of the ratio of specific heats. Measurements indicate that the mean mo¬ 
lecular weight is also nearly constant to an altitude of about 90 kilometers but 
decreases, in a fashion not well known, above this altitude. Actually then, the 
pressure, density and speed of sound measurements provide reasonably good 


Ionosphere: 

(continued) 

Mesosphere 

and 

Exosphere: 
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determinations of all atmospheric properties below 90 kilometers, the altitude 
range of most interest in this study. 

The U.S. Extension to the ICAO Standard Atmosphere, Reference 2.6, 
tabulates 17 important atmospheric properties to an altitude of about 90 kilo¬ 
meters and continues the tabulation to 300 kilometers for some of the 17. In 
the region from sea level to about 500, 000 feet, the variation of density can be 
described as roughly exponential; it decreases by a factor of approximately 
ten for each 54, 000-foot increase in altitude. Figure 2. 2 shows the density var¬ 
iation with altitude to 450, 000 feet as assumed by Reference 2. 6. Figure 2. 1 
shows the assumed piece-wise linear variation of temperature with altitude 
given in Reference 2.6. 


Density and Temperature Models: Of the physical properties of the at¬ 
mosphere, the density is the most important to the description of the aerody¬ 
namic forces acting on a re-entry vehicle. The lift, drag, and stability and 
control characteristics of the re-entry vehicle depend heavily upon dynamic 
pressure which is a function of density as well as the velocity of the relative 
wind. Similarly, the calculation of Mach number requires an accurately devel¬ 
oped model of atmospheric temperature. Each of thei*e physical properties is 
modeled as a function of geometric altitude to 300, 000 feet and the models are 
described in the following paragraphs. 


The basic source of data for both density and temperature is the U.S. 
Extension to the ICAO Standard Atmosphere (1958), Reference 2.6. This work 
tabulates atmospheric properties to 300, 000 geometric meters and correspond¬ 
ingly to 1,000, OOP geometric feet. There are other sources of such data; but 
for altitudes below about 300, 000 feet, there is very little difference among them. 
The density model is essentially based upon Table IV B, page 199; and the temp¬ 
erature model is based entirely upon Table IV A, page 185, Reference 2.6. 


Since above an altitude of about 300, 000 feet there is no meaning to the 
speed of sound and there is little aerodynamic heating for the slowly decelerating 
vehicle, the temperature model is developed accurately to an altitude of 300, 000 
feet and no higher. The speed of sound is needed to evaluate Mach number which, 
in turn, is needed to evaluate the aerodynamic coefficients; however, the aero¬ 
dynamic effects are so small at altitudes above 300,000 feet that accuracy in the 
temperature model is not important there. 


Figure 2. 1 shows the variation of kinetic temperature, 'T , in degrees 
Rankine with the altitude, ^ , in geometric feet. The speed of sound is defined 
as: 


Cl 



(2.9) 


in which, to an altitude of nearly 300, 000 feet, the following values hold: 

=1.4, the ratio of specific heats of air for 
altitudes 

/V = 28. 966, the mean molecular weight 


s 
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- 49719.64 ft-lb / *R slug, the universal 
gas constant 

~T = temperature in degrees Rankine 

This reduces to the expression programmed with the equations of motion: 

6L - 44./ {7* (2. 10) 

Because the temperature function is made up of straight line segments only, it 
is programmed in the IBM 704 at just the break points, and the IBM 704 linearly 
interpolates between these values which are tabulated below. 


h~ 


T~ 


Geometric 

Feet 


Degrees 

Rankine 


0 

518.688 

36,500 

389.988 

82,000 

389.988 

156,000 

508.788 

175,000 

508.788 

250,000 

354.35 


For altitudes above 250,000 feet, the constant value T = 354.35 *R is used. 

Since in this study the primary interest in atmospheric density lies in the 
region below 300,000 feet, special effort has been applied toward obtaining a good 
empirical fit of the standard atmosphere density data of Reference 2.6 below 
300,000 feet. The data is exponential in character; hence the empirical model is 
composed of three exponential segments, two of constant slope (on a plot of A 
versus log /> ) separated by a segment whose slope is a function of altitude. This 
approach allows the use of a fixed coefficient, the value of density at sea level, 
for the exponential function. 

The empirical model is plotted in Figure 2.2 along with a plot of the stan¬ 
dard atmosphere data. It can be seen that the model of density used in this pro¬ 
gram fits the standard atmosphere very well at altitudes below 300,000 feet. The 
plot is, however, extended to higher altitudes to show that the model is still good 
up to 450, 000 feet. The equations which is programmed and has been used to ob¬ 
tain the solutions of Sections V and VI is: 

eXyC?./_ 36 _ \ ( 2 - l 2 ) 

3/COO-z) 

0 £ h < 25000 ft 
25000 ft < h i 65000 ft 
65000 ft < h 

/O = 0.002378 slugs/ft 3 

The value for^e is the one that has been the standard for years and does not ex¬ 
actly agree with that of Reference 2. 6 (given as .0023769). Also, at h = 25,000 ft 


/* 


. (0. 


2 h ■ 
8000 
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and h = 65,000 feet, the empirical model has slope discontinuities, but the curve 
itself is continuous throughout. At the time this density model was devised, it 
was thought that an exponential having a slope which was a function of altitude 
would be the best way to fit the curvature of the standard atmosphere data. It is 
felt now, however, that three exponential segments, each having a different con¬ 
stant slope and a different coefficient, would be equally good. 

The atmospheric pressure is required for the description of rocket thrust 
in the atmosphere as in Appendix B. Pressure, P, may be obtained from the 
equation of state for a perfect gas: 

P= 3090.24 ^T lb/ft 2 (2.13) 

in which density, p , and temperature, T, are given by the models described in 
this section. Above 300, 000 feet the temperature is not significant to the speci¬ 
fication of the forces on the vehicle, so the temperature model used in this study 
(Figure 2. 1), as a matter of convenience, neglects any variation of T above 
300,000 feet. 

2.3.3 Winds and Gusts 

It is conceivable that a wind, such as a jet stream, could produce a sub¬ 
stantial change in a vehicle's path. It is conceivable also that turbulence, or gusts, 
could produce severe short period vehicle motions which would make piloting and 
control difficult. Each of these possibilities should be considered important to the 
over-all description of re-entry vehicle flight for a simulator. 

The motions of the atmosphere, called circulations, are described in Chap¬ 
ter VII of Reference 2.7. The primary circulation is that which exchanges air be¬ 
tween high and low latitudes. It is, by nature, a thermal circulation composed of 
three meridian cells of flow which is heavily modified by the earth's rotation and 
its radiative heat loss. This atmospheric circulation accounts for the general 
westerly surface winds in the middle latitudes and the easterly surface winds in 
the polar and equatorial zones. In the upper troposphere the wind direction is 
variable with season, being westerly in winter and easterly in summer (Reference 
2. 5). Differences in temperature between the oceans and the land masses, in cer¬ 
tain regions, cause circulations to and from summer-heated and winter-cooled 
continents. These are called monsoon circulations. Similar to these circulations 
are the minor circulations which result from a variety of local situations and are 
diurnal in character. The circulations described in this paragraph are always 
present because of the earth's surface features and the influx of solar energy, but 
they are often completely masked by large perturbations in the air flow associated, 
at the surface, with cyclones and anticyclones which produce daily weather changes. 

The westerly winds of the upper troposphere, mentioned in the preceding 
paragraph, prevail at all latitudes and are often concentrated in one or more jet 
streams. A jet stream is a pronounced concentration of west-to-east flow in a 
narrow band centered just below the tropopause in the middle latitudes around the 
earth. It appears as a fast-moving meandering stream whose position varies . 
considerably. The meandering can be associated with the cyclones and anticyclones 
of the air near the earth's surface. Of the two jet streams which are ordinarily 
present in the upper troposphere (two in each hemisphere),, the most northerly one 
(the one found near the border of the United States and Canada in the summer) is the 
stronger. 
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At much higher altitudes the data on wind velocities is very sketchy be¬ 
cause of the measurement techniques required. Basically, the summer winds in 
the region from 100,000 to 200,000 feet are easterly, averaging about 100 feet 
per second; and from 200, 000 to 300, 000 feet the winds gradually become west¬ 
erly, approaching a peak velocity slightly above 300,000 feet. In winter the winds 
are westerly, peaking at about 200 feet per second at an altitude of 200, 000 feet. 
Chapter VII of Reference 2.7 and Chapter V of Reference 2. 5 contain more in¬ 
formation about winds at these high altitudes. 

The description of the general winds presented in the preceding three 
paragraphs is simply an averaged condition based upon as much data as is avail¬ 
able. In the atmosphere at the earth's surface, and especially throughout the tro¬ 
posphere, the daily weather variations produce great changes in the wind velocity 
and direction at a given location. Thus, a general wind profile, as developed for 
this program is only a long-term average of the mean wind measured on many 
different occasions. 

Near the earth's surface at any particular time, the wind can be consid¬ 
ered a combination of a vector mean wind and the superimposed turbulence (gusts). 
To simplify matters, the vector mean wind is assumed constant over the field of 
interest, and gust velocity amplitude probability distribution is assumed Gaussian. 
This latter basic assumption allows the statistical description of turbulence. Stud¬ 
ies of data have shown that, except for very low altitudes, the assumption is ade¬ 
quate. In addition, indications are that turbulence takes the form of individual 
patches or fields throughout which it is reasonable to assume that spatial homogen¬ 
eity* and isotropy** hold. Caution must be exercised when using approximations 
of homogeneity and isotropy. For example, isotropy is a poor approximation for 
certain kinds of clear air turbulence and for turbulence under 1000 feet altitude. 

One patch or field of turbulence is different from another if the root-mean- 
square velocities and the velocity spectra are different. Some factors which affect 
turbulence intensities (velocities) and spectra are: 

1. Weather 

2. Terrain features 

3. Altitude above terrain 

4. Flight path direction relative to terrain 

5. Magnitude of the mean wind 

6. Direction of the mean wind relative to terrain 

7. Temperature gradient. 

This listing is not meant to be in order of importance of the factors as present 
knowledge of atmospheric turbulence is insufficient for such a classification. Ref¬ 
erences 2.8 and 2.9 give more details. 


* Spatial homogeneity means that the statistical properties (root-mean-square 
gust velocity and gust velocity spextrum) are the same at each point in die field; 
that is, the coordinate axes may be translated without effect. Temporal homo¬ 
geneity (stationarity) may be assumed when the statistical properties of the field 
do not change with time. 

** Isotropy means that the statistical properties are independent of angular or¬ 
ientation; that is, the coordinate axes may be rotated without effect. 
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There is considerable gust data for altitudes in the middle regions of 
troposphere, but for very low altitudes below a few hundred feet and for high 
altitudes near and above the jet stream, data is scarce. The most severe condi¬ 
tions are probably met in low altitude approach and landing, and though data in 
this region is meager, what there is should suffice for the simulation problem. 
Turbulence associated with jet stream shears, though not well known except that 
its intensity is high, is probably much less a problem because of the higher air¬ 
speeds and the less exacting pilot tasks at jet stream altitudes. 

Wind and Gust Models: The winds around the earth are variable in di- 
rection and velocity, dependent upon location, and dependent upon seasonal and 
diurnal effects. It is, therefore, difficult to do much more than develop a model 
based upon large-scale general winds though this may be somewhat unrealistic. 
The model presented below is developed so that it is a function of latitude and al¬ 
titude only. The wind direction is always the same; the function gives magnitude 
only. The bas^c source of wind data is Reference 2. 5 which presents data in tab¬ 
ular and graphical form. A gust model has not been programmed; however, a 
usable model is described. 


For the wind model, attention was centered upon mean winter winds in 
the 30 to 45 degree north latitude region. In this region the dependence of wind 
direction and velocity upon longitude is needlessly complex for simulator purposes. 
Hence, the directional data was ignored and the magnitude data was graphically 
averaged and used as if it were completely zonal; that is, in a direction parallel 
to lines of latitude. Data is presented for altitudes to 100,000 feet in the Vector 
Mean Wind Table of Reference 2. 5. For altitudes from 100,000 to 300, 000 feet, 
available data is much less complete and is given (Reference 2.5) in zonal form 
only. Here again, much averaging was done to put the data in useful form for this 
study. The result of the graphical averaging is the general winter wind profile 
shown as the solid line in Figure 2.3. The profile gives the wind velocity magni¬ 
tude with respect to the local surface of the earth as a function of altitude, h, for 
the 30 to 45 degree latitude zone. The wind is assumed to be always westerly 
(from the west) in both north and south latitudes. The profile accounts for the jet 
stream, and the general geostrophic winds at higher altitudes. 


An empirical fit for the wind profile is also plotted in Figure 2. 3 as the 
dashed line. This equation is: 

_ 17,280 

144 + (35 - .001 h) 4 1500 H 


V 4 a 30 + • 


300 » 00 °-J|.in^| (2 - 

) + (20 - . 0001 hr 


14) 


where is the westerly component of the wind velocity (ft/sec) 

with respect to the local earth's surface. 


The wind velocity matrix [W] has northward and upward components also, which 
in this study are assumed identically zero. The absolute value of the sine of 
twice the latitude angle produces the variation of wind velocity with latitude. The 
magnitude of the wind velocity is zero at the equator, reaches a maximum at 45* 
north or south latitudes, and is zero again at the poles. This variation is in rough 
agreement with the measured variation of zonal wind magnitudes at various latitudes. 


No gust velocity model was programmed in this study, but such a model can 
be devised by suitably modifying the output of a digital random function generating 
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subroutine using some digital filtering technique to obtain the correct gust fre- 

quency spectrum. The filtered output is then scaled to the desired level and 

used directly as yu x ,, and. 
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FIGURE 2. 2 ATMOSPHERIC DENSITY MODELS 
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FIGURE 2. 3 ATMOSPHERIC WIND PROFILE MODELS 
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SECTION III 

EQUATIONS OF MOTION 

3.1 INTRODUCTION 

The general purpose of this, study is to develop equations of motion suit¬ 
able for simulating the re-entry and landing problem for a manned hypervelocity 
boost glider, and particularly for application to flight simulators using digital 
computers. The equations of motion as derived for use in this study, as distin¬ 
guished from those suitable for a simulator, accent the new dimensions assumed 
in dealing with space flight and hypervelocity vehicles; and their purpose is to 
determine how these new dimensions must be handled. As such, these equations 
embody many new concepts not normally found in flight simulator equations. On 
the other hand, many characteristics commonly included in simulator equations 
will not be found in the equations used in this study. It is not the purpose of this 
study to investigate matters which are already well understood, but rather to in¬ 
vestigate the new facets being introduced. The complete equations of motion des¬ 
cribing the flight of boost glide space vehicles, including both the attitude and tra¬ 
jectory problems, are exceedingly complex. Thus, this study is essentially lim¬ 
ited to piloted flight - the over-all re-entry and landing problem - and emphasis 
is placed on determining what simplifications can be made in the complete equa¬ 
tions with the piloting task in mind. 

The equations of motion as derived for this study include the following 
characteristics and limitations: 

1. Six rigid-body degrees of freedom for the over-all vehicle, 

2. A spherical earth rotating with constant angular velocity, 

3. Flight limited to sub-orbital speeds and below 300, 000 feet 
of altitude, 

4. No restrictions on the vehicle's attitude, 

5. No restrictions on the vehicle's trajectory within (3) except 
that flights directly over the North and South Poles are 
excluded, 

6. Non-varying maps (i. e. , limited power), 

7. Aerodynamic and reaction controls, 

8. Non-linear aerodynamics, 

9. Undisturbed atmosphere rotates with the earth, 

10. Winds defined in terms of northward, eastward and upward 
components as functions of altitude, latitude and longitude, 

11. Gusts or turbulence defined in terms of the vehicle - lateral, 
normal, and head-on components, 

12. Gravitational force field for a spherical earth, and 

13. Air density and temperature variable with altitude. 

None of the restrictions listed will compromise the objectives of this study, 
but the restriction to a non-varying mass (6), and to a spherical earth (2) and (12), 
might seriously hamper extension of the equations to the full flight of sub- or near- 
orbital vehicles. Accordingly, methods for including the rapidly varying mass 
and large thrust forces associated with boosted flight are developed in Appendix B, 
and methods for including earth oblateness and surface aberrations are developed 
in Appendix C. 

The restrictions placed on the flight range are compatible with atmospheric 
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X 


flight, but preclude flight where the aerodynamic forces are negligible. For 
flight outside the atmosphere, many effects become important (e.g., gravita¬ 
tional aberrations, radiation pressures, magnetic fields) which are missing in 
the equations used in this study, and which even the smallest aerodynamic forces 
and moments would completely mask. The restriction to non-polar flight is an 
artificial one produced by using latitude as a variable. If polar flights are re¬ 
quired in any particular application, one of two methods developed in this study 
can be used to remove this restriction: (1) the method used to handle a 90* side¬ 
slip, which involves using direction cosines instead of angles in some relations, 
or (2) the method used to handle vehicle attitude, which involves switching to an 
alternate axis system. 

The general development of the equations begins by setting forth basic 
principles and concepts, which are formulated most conveniently in vector nota¬ 
tion. Next follows the selection of the axis systems, coordinates, and transfor¬ 
mations which are necessary to develop the equations in detail. Finally the de¬ 
tailed equations are developed in matrix form. The extensive use of rhatrix al¬ 
gebra is dictated by its concise form and the ease with which matrix equations 
can be incorporated into a digital computer program, the method selected for 
solving the equations in this study. Since the algebraic form of the equations is 
of some interest, the matrix equations have been expanded and are presented in 
Appendix A. 

3. 2 VECTOR EQUATIONS - BASIC PRINCIPLES AND CONCEPTS 

This section is concerned with the application of the laws of mechanics to 
the dynamical description of the vehicle, and its relative motion with respect to 
the earth and the atmosphere. The fundamental assumption is that the vehicle 
may be treated as a rigid body. 


3.2. 1 Equations of Motion 


The derivation starts with the assumption of Newton's Laws of Motion for 
a particle and their extension to a system of particles. The inertial frame is 
taken as fixed with respect to the "fixed stars", and secondary (moving) reference 
frames are defined by the earth and the vehicle. Coordinate systems and position 
vectors for these reference frames are shown in Figure 3. 1, and it should be 
noted that X 0 . / . and g are parallel to , Y r > and Z z , respectively. 


The basic vector equations of motion for a body or system of mass par¬ 
ticles, where the moments of forces and moments of momentum are restricted 
specifically to be taken about the c. g. , are written following Reference 3. 1 as 


m-i 

m • * 


where the derivatives are taken with respect to inertial space as denoted by the 
"I" subscript. The forces exerted by the sun, moon and stars - called the as¬ 
tronomical forces - are separated from those exerted by the earth as 
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Tj = F + F astronomical 
Mj = M + M agtronorri i ca i 


Treating the body as a constant mass, and writing the force equation in 
terms of the position vectors, we obtain 


/ cffii \ •• •• •• — 

y tLf /f = + '** ^astronomical 


(3.3) 


Consider the force on the body (vehicle) due to the sun, the largest astronomical 
force. Applying the inverse square law of gravitation yields 


(U 


F sun (on earth) 
m (of earth) 


i £S - 


F sun (° n v^icle) />j\* - 
m (of vehicle) ^ S£ J 


where is a unit vector in the earth to sun direction 

Jgg is a unit vector in the body to sun direction 
Y* tg is the distance from sun to body 
is the distance from sun to earth 
F^ w is the force exerted by the sun. 

Since (Y\ t /r\ e f-l and =/^, then 

W^) r = F sun < on vehicle) = F sun 

Extending the results to all other astronomical bodies yields: 




astronomical 


By definition, the origin of the body axes are the c.g. , so that 
Mr = M since Mastronomical = 0 


(3.4) 

(3.5) 


Now, subtracting Equation (3.4) from (3. 3) and using (3. 5) in (3.2) yields the 
following: 




(3.6) 



* The equality is exact for a spherical mass distribution. 


(3.7) 
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Thus, for motions near the earth, we can treat the earth's center as an inertial 
reference provided the gravitational forces of all celestial bodies other than the 
earth are excluded. In this manner, X 0 , Y , and z* can be considered the 
inertial frame. * 6 13 

The equations of motion can be written in more explicit form using body- 
reference quantities as follows. 


V 

(3.8) 

W = I to 

(3.9) 

where I is the moment of inertia operator (Reference 3.1)** and the inertial 
frame is the pseudo-inertial earth-centered reference frame. The general equa¬ 
tions of motion for the vehicle when treated as a rigid body follow. The equations 
are extended to deal with variable mass in Appendix B. 

+ to x ~?rv V - F 

(3.10) 

I £3 +■ IS k X to = /V 

(3.11) 

vdri - % 

(3.12) 

CO = f £> dT+ tO c 

(3.13) 


(No subscript is needed on u> since (co) T - (u))& ). 

The integrals for V and LO must be interpreted with care, since general vector 
integration is not definable (Reference 3.2). These integrals imply integration 
of the components of (V) e and O) in the body axes. 

3.2.2 Kinematics Relative to the Earth 

It is necessary to describe the trajectory of the vehicle and the rotation 
of the vehicle relative to the earth, that is, with respect to the X^-, Y g , 
frame which rotates with, and is fixed in, the earth. The relative linear and 
angular velocities of the vehicle may be expressed as: 

V e = V- (3.14) 

U t - u - JL (3. 15) 


* This does not imply that the earth is travelling in a straight line and is a 
true inertial reference. The approximation is much better than that. It im¬ 
plies that the acceleration of the earth and the vehicle due to astronomical 
forces are equal, that the acceleration of the vehicle relative to the earth 
arises entirely from earth forces. 

** T is a dyad or second order tensor and Equation (3.9) would be written as 
X- T-ss in the notation of Reference 3.2. 
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where the velocity due to earth rotation is 
The position of the vehicle is obtained from 

ft.); - Z 

so that 



(3.16) 

(3.17) 


(3.18) 


The integral for implies integration of the components of Vg in earth axes, 
and 9* ee is a vector fixed in the earth axes which denotes the initial position 
of the vShicle at t = 0. 


3.2.3 Kinematics Relative to the Air 

In order to determine the aerodynamic forces on the vehicle, it is neces¬ 
sary to describe the motions of the vehicle relative to the air. The undisturbed 
air is assumed to rotate with the earth, and two types of disturbance velocities 
are superimposed on the general rotating atmosphere: 

1. general air mass motions or winds, W , which are defined 
with respect to the earth, and 

2. rapid air mass motions or gusts, Z , which are defined 
with respect to the airplane. 

The reason for this duality is that the accelerations and the velocities of the air 
relative to the vehicle must be specified in order to determine the aerodynamic 
forces. This is an extremely complex problem where the winds are specified 
relative to the earth, as it requires determining nine velocity gradients with re¬ 
spect to the airplane (i. e. , 9 W x /dX , 9Wy/dX , .-‘.dW^/dy , •••» etc.) and 
their time derivatives. Where the winds are expressed in terms of empirical 
geophysical data, this represents a fantastically complicated computing problem. 
In order to avoid this situation, several assumptions are made: 

1. The air motions can be separated into low frequency motions 
(winds, W ) and high frequency motions (gusts and turbulence,^ ). 

2. The gradients of the low frequency winds are small enough to 
be neglected. 

3. The high frequency air motions, gusts and turbulence, can be 
specified directly in terms of components in the body reference 
frame. 

These three assumptions rest on relatively firm ground. A fourth assumption, 
not so firmly grounded, that one would like to make is: 

4. The relative accelerations of the low frequency winds, and also 
those due to the rotating atmosphere, may be neglected in deter¬ 
mining the accelerations of the air relative to the vehicle. 

This fourth assumption will be investigated in the course of this study to see 
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if it really is valid. 

The velocity of the vehicle relative to the air (measured at the c.g.) and 
its time rate of change in the body axes frame are given as: 

V A = V - V er - W -JL (3. 19) 

‘ &) ES -(ji) B + (3.20) 

in which u>g. e$ is the angular velocity of the body axes frame relative to the earth- 
surface frame. The angular velocity of the airplane with respect to the air (mea¬ 
sured at the c. g.) would be given, in the absence of any assumptions, as 

= uj - cS w - -JL 

However, assumption (2) implies that = 0. Furthermore, the present knowl¬ 

edge of atmospheric turbulence does not provide any data on £5*,, so that perforce, 
- 0. Thus we are led to the equation 

a <2 e = & - J1 (3.21) 


3.2.4 Forces and Moments 

The forces and moments acting on the vehicle can be grouped functionally. 


Fa. 

\ 

aerodynamic forces and moments 

5 


gravitational forces and moments 

F r 


thrust and engine forces and moments 
(power), and reaction controls 

af 

h/7 

other forces and moments. 


The aerodynamic forces are generally functions of air density, air velocity, and 
aerodynamic coefficients. The gravitational forces are generally functions of the 
position of the vehicle with respect to the earth. and are forces and mo¬ 
ments expressed directly in dimensional units (lb and ft-lb) and are functions of 
some arbitrary control. &V and are reserved for additional terms which may 
arise, say from a booster or rotating machinery. 

3.2.5 Measured Accelerations 

The accelerations that would be felt by the pilot or measured by acceler¬ 
ometers (three components) placed at the c.g. of the aircraft are given by the 
external forces (i. e. , they do not include the gravitational force) as follows. 

■W = (F-r,)/** * (£ 

Now, if the accelerations are measured at some place other than the c.g. of the 
aircraft (as would usually be the case), the the effects of rotation must be included. 
If the instruments are located with respect to the c.g. by the position vector ^ , 
then the acceleration of is given as: 

(r^ = (r t -) e + u> x y\ + z & x Cfy) a + eS x(Z5x r> ) 
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For a rigid body (no flexibility) 

OH - (H ‘o 

Thus ^ , the measured acceleration, is given as 

■ (F - Fq + to x + u> x (coxr\) 

3.3 AXIS SYSTEMS, COORDINATES, AND TRANSFORMATIONS 

Matrix notation is fundamental to the description of coordinate systems 
and their transformations, and hence, is used extensively in this development. 

The assumption of a spherical earth is inherent in the development, but the ap¬ 
proach can be extended to include oblateness as outlined in Appendix C. 

3. 3. 1 Earth-Surface Axes 

Spherical coordinates have been selected for defining the position of the 
vehicle with respect to the earth. Such coordinates have an obvious disadvantage; 
a singularity exists at the poles and precludes polar flights. On the other hand, 
these are the coordinates commonly used in both terrestial navigation and in ce¬ 
lestial mechanics; and, as will be discussed in Section 3.3.2, they are the means 
for removing an even more emharassing singularity. Accordingly, "earth-surface 
axes" ( X £S > ’ and 2 t% ) are defined as shown in Figure 3.2. The X £ $ and y es 

axes are parallel to the earth's surface (directed north and east, respectively), 
and the « axis points toward the earth's center; but the origin of the earth- 
surface axes is located at the vehicle's c.g. and not on the earth's surface. Since 
the earth is considered spherical, the orientation of these axes is described by 
the longitude ( A ) and the latitude (^ ) of the vehicle, and their position is des¬ 
cribed by these same quantities and the magnitude of the radius vector < *)n )• 

The position could be described alternatively by components in the earth axes 
and^, ). 

Convenient relative velocity components are defined by the earth-surface 
axes: the vehicle's velocity north ( V„ ) and east ( V*. ), and its rate of climb 
( ]/ ). The relative tangential or horizontal velocity vector ( \? r ) is described 

in magnitude ( V_ ) by the vector sum of ty and , and in direction by the 

flight-pathjieading angle ( Tjr T )• The flight-path angle ( 7* ) is the angle between 
V r and V e . 

3.3.2 Orientation of Body Axes 

Conventional body axes ( % , y , 2 ) are selected to represent the vehicle. 
The method of orienting these axes in conventional airplane dynamics is well es¬ 
tablished (e.g. , Reference 3. 3): three sequential rotations or Euler angles ( fr , 

6 , and ) are used to rotate from inertial axes to body axes. This same con¬ 
ventional scheme has been used by Patha in Reference 3.4 to handle the larger 
problem associated with hypervelocity vehicles, but it leads to serious difficulty. 
A singularity occurs whenever the X -axis is parallel to the polar axis, a likely 
situation in any re-entry and landing trajectory. Patha solves this problem by 
re-orienting the fixed reference frame each time a singularity is approached; 
but in so doing, the orientation angles lose their physical significance. The ap¬ 
proach taken here is to introduce additional orientation angles (five in all), thus 
strategically locating the singularities where they will do the least harm. In ad¬ 
dition, an alternate set of orientation angles is resorted to whenever the 
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singularities are approached but once passed, the original angles are restored. 
This scheme was suggested by Patha's re-orientation scheme, but has the addi¬ 
tional advantage of retaining the physical significance of the angles. 

The usual definition of the orientation angles ( p* , & , and 0 ) is made 
as shown in Figure 3.3a, except the reference frame is the earth-surface axes 
rather than the earth or inertial axes. Thus the usual significance of these angles 
is still retained: pr is the heading from north, Q is elevation from the local 
horizontal, and 0 is the bank angle from the local horizontal. 

The singularities occur for a vertical attitude ( O = *90*). Since Q - 90* 
is a possible re-entry attitude, it is deemed necessary to provide for this case, 
and also concurrently the case of 0 = -90* which might occur in a spin. Ac¬ 
cordingly, a second set of "special orientation angles" is defined ( » 6> s . 

and 0 S ) as shown in Figures 3. 3b and 3. 3c, which is used only when Q —► -90* 
(nose down) or Q —> +90* (nose up). These angles can be thought of as being the 
same angles as the normal ones, except the reference position of the vehicle is 
with the nose pointed straight down. Alternatively, these angles can be thought 
of as providing the azimuth of the pilot's head ( ^ ), the elevation of the pilot's 
head ( ), and the wing tilt angle ( 0^ ). In terms of nose position, A 

gives the upward tilt from straight down ( 0 near -90*) or the downward tilt 
from straight up ( & near +90*); 0 ^ gives the leftward tilt from straight down, 
and 0 S - 180* gives the leftward tilt from straight up. These are relatively sig¬ 
nificant quantities when the attitude of the vehicle is vertical (nose down or nose 
up). But it is perhaps more important that this particular selection of the alter¬ 
nate "special orientation angles" is most easily incorporated into the equations of 
motion. 


The method for switching from one set of axes to the other employs a cri¬ 
teria based on & . When Q reaches some large negative or positive value, say 
|@| = 70*, the coordinates are switched from 0,0, and pr to 0^ , 6> s , and 

. Then when 9 again reduces, say to |0| = 50*, the coordinates are switched 

back. Different switching values of |©| are used so that small oscillations will 
not cause continual switching of coordinates. 

The angles pr , Q , and 0 (alternatively ^ , and 0 % ) are only 

part of the general problem, orienting the body axes with respect to the earth- 
surface axes. It still remains to orient these latter axes with respect to the 
earth axes and the inertial axes. The earth-surface axes and the earth axes are 
related by the longitude ( X ) and the latitude ( ■p ); the earth axes and the iner¬ 
tial axes, by the angle -At . The angle -At is only of interest if celestial navi¬ 
gation is specifically included. The complete orientation problem (for normal pr , 

Q , and 0 ) is depicted in Figure 3.4. The sequence of rotations (orthogonal 
transformations) going from earth to body axes is as follows: 

1. a rotation X about to X t > Y > £, 

2. a rotation "0 about - to X* > Y z > 

3. a -90° rotation about Y f to 

4. a rotation about ^ <s to X, , 

5. a rotation 9 about ^ to X z , 

6. a rotation 0 about X t to X , , J. 
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To complete the picture, the sequence for the special orientation angles 
( , © s , 0 g ) is unchanged through rotation (3), and continues as: 

4. a rotation ^ about to ^ 

5. a rotation © s about y^ to X* , y^ • 

6. a rotation 0 g about X* to X 3 , 

7. a - 90 ° rotation about y^ to X, y , 

3.3.3 Wind, Stability, and Body Axes 

The wind axes and stability axes are defined in the conventional manner, 
as shown in Figure 3.5^ The only unusual characteristic is that these axes are 
defined by the vector , as distinct from the usual vector V . The body axes 

are located with origin at the c.g. and % and are in the (assumed) plane of 
symmetry; but the "incidence" of X is essentially arbitrary, and may be set 
at any convenient value. The sequence of rotations (orthogonal transformations) 
relating the axes, going from wind to stability to body axes, is as follows: 

1. a rotation/3 about <^^to •fa 

2. a rotation about y^ to X , y , 

3.3.4 Coordinate Transformations 

It is necessary to define coordinate transformations between the var¬ 
ious axis systems. These are defined in terms of the direction cosine ma¬ 
trices, and differentials are used to avoid including origin shifts in the basic 
definitions. 


Earth axes to earth-surface axes: 




*Y t 


Earth-surface axes to body axes: 

T d* 

dy 
** 


Earth axes to body axes: 


and 

Wind axes to body axes: 



* M 

' [".-«] 
- M 


d *« 
d Ui 

"f. 

?•«*] 



(3.22) 

(3.23) 

(3.24) 

(3.25) 

(3.26) 
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The foregoing transformations are orthogonal, each matrix element 
being the cosine of the angle between the corresponding axes. These can be 
handled directly in terms of direction cosines, but the approach taken here is 
to express the direction cosines in terms of specific angles - Euler angles or 
orientation angles - which have some general physical significance. The 
orientation angles are defined in the previous section, and it remains to define 
the transformations in terms of these angles: 

M - WWW 
['•».«] -WWW 
[<•».*] *HM « 3 «» 

And the transformations representing the individual rotations (single angle 
transformations) are listed as follows: 


(3.27) 

(3.28) 


[*] 

M 

[•] 


cotX 
- tin X 
O 



cos e 

° 

0 

[ cos JB 

r 


«/<r 1 O 

cot X o 

O / 


o 

'1 

o 


o 

0 J 

o 

-tin ©1 

1 

0 

o 

cos ej 


- fin J& 0 
cos JB O 
O I 


M 

W 

(•) 

W 


|~eos i> o sin 1>~| 

° ' ° 

o cot 

E co* ir sm ir o “1 
•tmir cot if o I 
O O ; J 

p 0 ° ] 

O CO* # sin 0 I 
[_0 -tin 0 cos 0J 

[ cose O 

0/0 
tin «. o cose 


The definitions are not complete without specifying the range for each of the 
angles: 


-180* 4 X, £ 180" 
0 4 p 360* 
-180* 4. fi * 180" 


-90* * V £ 90* 

-90" ^04 90"; 
-90* * u 4 90" 


-180* * 0 4k 180? 


In addition, the transformations for the "special orientation angles" 
must be defined. The symbol [i. 9 J is used to describe the matrix when 

this matrix is defined in terms of the special orientation angles ( ^ , O, , ). 
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N ■ K„) 

M 3 MWWW 


(3.30) 

(3.31) 


The functional form of the w w and matrices is identical to that of 

the M . [©] , and[^] matrices (the trigonometric functions of cor responding 
elements), and the range of ^ , © s , and 0^ are the same as those for "iff , 

Q , and 0 , respectively. 

Finally, it should be noted that the stability axes are the intermediate 
axes between the wind axes and the body axes. 


(3.32) 


‘ dx' 


~dx % ' 



■M 

[A,l 

dtf 

• H 


or 

<*** 


_ d f_ 


/V 


d % 




3.4 


MATRIX EQUATIONS 


The basic principles and concepts have been developed, using primarily 
vector notation. Axis systems have been selected for describing the various 
physical quantities, and methods for transforming quantities from one axis sys¬ 
tem to the next have been defined. The next step is to write the detailed equa¬ 
tions describing the kinematic and dynamic behavior of the vehicle. Matrix al¬ 
gebra is used almost exclusively, in order to keep the development tractable, 
and because matrix equations are easily adapted to digital computers. 

3.4. 1 Equations of Motion 

The equations of motion are written in body axes. The linear and angular 
momenta (from Equations 3.8 and 3.9) are written as 

(3.33) 

(3.34) 


M * «M 
M - MM 


M 


. M 


• M 


[a;?] 


The form of JJ J implies that is a plane of symmetry. 

The notation generally adopted is that the primary symbol refers to a 
specific type of quantity (e.g. , V for the velocity of the body), and the subscript 
indicates the reference frame involved (e.g. , Vg for the velocity of the body with 
respect to the earth axes) - no subscript usually implying an inertial reference. 
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I 


» 


The equations of motion are obtained (from Equations 3. 10 and 3. 11) as 

HI [a] - w [v]*m [ctfx] [v] « [r] (3.35) 

MM * WMM * M < 3 . 36 > 


where 



V 



M* 


. H 1 


. M- 

M v 


l a t\ 


L'»J 




and Ia)X 
3.5). 1 ' 


an operational form fo* - 

~ O 

[o?x] a 1* 

f 


the vector cross product 2>X (see Reference 



Note that the dot (e.g. , in [v] ) indicates a scalar differentiation of the individ- 
ual elements, not a vector differentiation. 


The above equations (3.35 and 3. 36) can be expressed in a form suitable 
for integration as follows. 

lv).i[r]-M[v] 


W-M'tN-MMM) 


(3.37) 
(3. 38) 


where 





° W z ] 

i/i, o , r-r,r -j* 
o z t /i J 


The velocities are integrated in the body axis system (following Equations 3. 12 
and 3. 13) giving 


[v] • ^ [v] rf7 ♦ [ V,1 (3.39) 

[«>] - f [v] it* [«»] 


3.4.2 Trajectory Kinematics 

The motions of the vehicle relative to the earth are described primarily 
in terms of quantities measured in the earth-surface axes (see Figure 3.2). The 


a 
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relative velocity of the vehicle is written (from Equation 3. 14) as 

M-KJM-K.) 


where 

"V4" 

[vj» 

"o 

kl - 

Vv 

V*r 


-v; 


0 


and the prime denotes the transpose (i.e. , for an orthogonal matrix the inverse 

is equal to the transpose; e.g. _| »|f._ _J We find, from Equa 

tion 3. 16, that L * S ‘* J 1 * “ J L 

V *r*nr £B *°* V (3.42) 


The magnitudes of the tangential velocity ( V r ) and the relative velocity ( Vg ), 
and the azimuth of the tangential velocity ( lU r ) and the elevation of the relative 
velocity or the climb angle ( y ) are given by 


V r */V 

i Vr>0 

(3.43) 

v e ■ /v/ . V/' 

i v e >o 

(3.44) 

% ±ir r *360* 

v r 

Quad. I II III IV 

Vy + - - + 

V A + + - - 

| (3.45) 

\A. ri 90* for- V r >0 

T * ; 

V 0 * T ±- 90“ for V r * 0 

(3.46) 

The position of the vehicle is described in terms of the latitude ( V ),Jthe longi¬ 
tude ( \ ), and the altitude ( k )• These are obtained by integrating Vg as in 
Equation 3. 18, but in spherical coordinates. 


rfT *[*«-<.] 

(3.47) 
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k = V r 

(3.48) 


55 »•« 

(3.49) 

. V, 

A =- cos V 

r *B 

(3.50) 

*•«•***. 

(3.51) 


3.4.3 Vehicle Orientation 

The angular position of the vehicle is described by the various orientation 
angles; specifically the orientation of the body axes with respect to the earth- 
surface axes is specified by the Euler angles , © , and 0 (or the special 
angles ^ , Q t , and 0 S ). These angles are obtained from the angular velocity 
relationships. 


"Normal" Euler Angles: The angular velocity of the body axes can be ex- 
pressed as the sum of the angular velocity of the earth axes ( 7i ),_the angular 
velocity of the earth-surface axes with respect to the earth axes (tOgg-t)' and 
the angular velocity of the body axes with respect to the earth-surface axes 

<"»-«>• _ _ 


Further subdivision can be achieved by introducing the angular velocity vectors 
representing the time rates of change of the orientation angles. 

z> e •[&+ 9 + ?) *($+i) + n ( 3 . 52 ) 

It is convenient to introduce the angular velocity of the vehicle relative to the 
earth ( 0) £ , as defined previously by Equation 3. 15) so that 

- & + $+++ X (3.53) 


Now this vector equations can be written in matrix form, transofmring 
ponents to body axes, as follows. 


M 


‘M 


MWiHO -vU'IM 


0 ) 


all com - 


(3.54) 




(3.54a) 


f 
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K1 - 


Equation 3. 54 can be solved for the rates of change of the Euler angles as 


[**-«]■ & -[*] H ■ M '£ 


w ■ mw 

W ■ 


<3.55> 




//cos & 
o 

sm e/cos 90/ 


o o 
/ o 


Finally, the desired Euler angles are obtained by integrating Equation 3.55. 

r*i 


[**-«] 




(3.56) 


The angular velocity of the body axes with respect to earth axes is obtained by 
writing Equation 3. 15 as 


n - m -i^m 


(3.57) 


"Special" Euler Angles: The "special" Euler angles are treated in exactly 
the same manner as the normal ones, with the following results: 




■NMN-WW 


(3.55a) 


N- W M'M' 
W= [«0 MMM 
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M £ ».] 


(3.56s) 


The functional forms of the "s" matrices are identical with those for the normal 
system, and the only change that is introduced in the equations is the insertion of 

the K matrix in the equations for . Of course, though the functional 

form of is different from ^ 5 j » these two matrices are actually equal 

by definition. 

The similarity between the equations and the matrices associated with the 
two sets of angles makes it very convenient and easy to switch from one set to the 
other; and it is for this very reason that 0s • , and are defined as they 

are. In a digital computer the same storage locations and the same operations 
can be used for the normal and the special sets. It is only necessary to make the 
one modification going from Equation 3.55 to 3.55s, and to compute the "initial 
conditions" required for Equation 3.56 or 3.56s. 

Interchange Between "Normal" and "Special" Angles: As mentioned in 
Section 3. 3 .1, the switching criteria are based on the value of & , and two values 
are used in the criteria. 


» 


0 - the value of 9 for which the change is made from the 
' 0,9 , ft set to the 0, , % set. 

Q- the value of 9 for which the change is made from the 
0 S , 0 9 , ffr set back to the 0 , 9 , ft set. 

In order to avoid continual switching due to small oscillations. 


\e M \ < \o,\ (3.58) 

The two matrices, and £/. s j , are written out in partial detail in order 

to clearly show the relationships used to calculate the starting values when 
switching occurs. 


M - 


COS & cos ft 


COS 9 sin it 
t-ZZ 


- sin 9 
sin 0 cos 0 
COS 0 COS 9 


M - 


[-COS 9 f COS ft. 


-Cos e. Sin ft. 


cos 0 t cos e % 
s/n 4 cos e, 
sin 9 


% 
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It should be noted that the elements in these two matrixes are always equal, 

L V 

When switching from the 0 , q , p* set to the . Pl set, 

the values of 0 t , , and ft are known and it is necessary to compute the 

corresponding values of 0$ > Q t , and . These values are computed by 

equating appropriate elements of the A and matrices. 


^ » sin ’'(cos 0, COS ©J ■ s/V'Zj, ; ° < < *°* L99 > ° 

O * >-9o\ to* L„<C 


(3.59) 





sin m/ 



f <?uap . I a nr ir\ 
{ L » ’ 1 


t>\ 


(3.61) 


When switching from the ~ n. , 0 % , ft set to the 0 , Q , set, 

the values of 0 , 6^ , and are known and it is necessary to compute the 

corresponding values of . B z . and P'z • The approach is the same as before. 


•i « •*» (-*>» J i 

< 90* /•* L i/a < O (3.62) 

O > S>z * ~to*0o* L $fs ^ O 



The switching time for a continuous solution of the equation, such as 
would be obtained on an analog computer, would be that time at which |0| =|©,| 
or \&\ = Kl occurred, depending on which set of angles was being used at the 
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time. For digital computation, however, it is more convenient to define a 
test value of Q , and then switch coordinates at the end of the first time inter¬ 
val for which the test shows switching is needed. Also, since 

= L t„ * -*»»>« 

it is more convenient to apply the test criteria to L, 9 and rather than © 

directly, as follows: 

(a) Switch from 0 , 0 , pr to 0 f , , ph when 

|^/*| ^ 

(b) Switch from ©L , 0. , 2t to 0 , 0 , 0r when 

The values of 0 /y and Q MT which have been selected in the present study are 
®, r = 70\ © jr =50*. 

3.4.4 Relative-Wind Kinematics 

The relative wind is defined by the inertial velocity of the vehicle and 
by the motions of the air due to earth rotation, winds and gusts. The specifi¬ 
cation of these quantities follows from the development in Section 3.2.3. 

The velocity of the air due to earth rotation, as measured at the origin 
of the earth-surface axes and in terms of components in these axes, is 



where 


l^ r = SLr rm Cot 1> 


The general air mass motions, or winds, are also measured with respect to 
the earth and are expressed in terms of components of earth-surface axes as 


M - 



(3.66) 


In the given order, these components represent northward, eastward, and up¬ 
ward (*W r ) air mass velocities. The rapid air mass motions, or gusts, are 
defined with respect to the vehicle and are, therefore, expressed in terms of 
components of the body axes as 


M 



(3.67) 
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The components represent, in order, the head-on, lateral, and normal gust 
velocities. * 

Using the three quantities defined above, the velocity of the vehicle with 
axel^a* t0 a * r * S °^ ta * nct * * rt>m Equation 3. 19 and with components in body 

[K] - [-* 

in which the matrix t-g. es transforms V er , and W from earth-surface axes 
to body axes. The airspeed ( V a ), the angle of attack (<X ), and sideslip ) 
are given by 


* M " { [ V„] * [*]} -|>] (3.68) 


K * r + w£' i |£ 


(3.69) 


„ [<?UAD. X U ju jt) 

* *<» —±- i { v- - - y- ) 

^ U — - -/ 


(3.70) 


o< ai s. 9o* /.* ur^ >o 

; Z.-90* £•* ur^<0 


K Atm o Ad /9 ± TO •* 


(3.71) 


Figure 3.5 shows the velocity vector, V?, , and indicates the reference directions 

for positive <* and /S . 


In Section IV, it is shown that aerodynamic representation of the vehicle 
requires, in addition to , « , and /9 , the variables Ot and ft which are 
used with the derivatives associated with non-stationary flow or lags (downwash 
lags). The solution for ai and £ is one of the stickiest algebraic problems in 
the whole development, and has been neglected by other authors (e.g. , Reference 
3.4). To solve for d and /& , the air velocity in wind axes 

W‘[s]“W N 

is differentiated and arranged to obtain: 
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§ 




< 

Vo.? 


Ar^ 

V x ot cos/3 



To find , Equation 3.20 is rewritten in slightly different form as 

(fog - &)g - (?„)' - {p) £s - Z e *-c ' & - We *■ &) (3.20a) 


When written in matrix form, Equation (3.20a) becomes 

H = W -[«•—«] {[(■„.«] [dUj ♦[*]) (3-73) 

-M [<w„ *] M-[>]*["« *] 

Certain matrices in this equation are specially synthesized 

[ Cu, *s-c)e S *] 

2 ' [ L *9-m ] [ ( 

and both are defined in the List of Matrices near the front of this report. Now 
and QL cot/3 can be expressed explicitly from Equation 3.72 in terms of 
, and tir^ as 

^ * { n/t Cat * *" Cot/9 sin/3 sin Of) (3.74) 


d cos/9 


(- u A sin « + cat «) 


(3.75) 


and the velocity components are obtained from Equation 3.73. The quantity 
QL cos /9 is obtained, instead of <& , in order that Equation 3.75 may still be 
evaluated at /3 - *90*. 


The angular velocity of the vehicle relative to the air is expressed from 
Equation 3.21 in body axes as 

[fcj * ["] “ H'M (3.76) 


M* 


The components of 
vehicle. 


**] 


are needed for the aerodynamic representation of the 


) 
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3 .4.5 Forces and Momenta 

Following the general description of vector forces and moments given 
in Section 3. .1.4, Iho forces and moments can be expressed in body axis com¬ 
ponents as 


M * K ^ [ F f] +[*] + [* F ] (3.77) 

M = W +[» 9 ] (3 - 78) 


For present purposes, we assume there is no thrust, no gravitational torque 
about the c.g., and hence set 

[/>] -o, [aAf]»o 

Appendices B and C provide methods for handling the forces and moments which 
by assumption have been neglected here, especially those due to thrust and 
gravity. 

Aerodynamic forces and moments are expressed in the usual coefficient 
form in body axes 


[*]■*/> (3 ' 79 > 

(3 - 80) 


where 

P 

s 



The matrices 
iables /V , * , 0 , 


= air density, a function of altitude given 
by Equation 2. 12 

= airspeed, given by Equation 3.69 
= the vehicle's wing area 


= the characteristic length matrix with reference 
span b , and reference chord sC 


= the aerodynamic force coefficient matrix in 
body axes 


= the aerodynamic moment coefficient matrix in 
body axes. 

and , in general, are complex functions of the var- 

<£ . fi , S e . $1 . Sp » '$*’**■ ancl others - There 
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is a peculiarity in the digital computer program, however, which precludes 
the dependency of upon a and ft . This arises because of computation 

order in the program: is required to compute [ft] which, in turn, is re¬ 
quired to compute [\JrJ Then [l/J is used to compute [li] and hence A 

and ft . Thus, though at and ft may be used to compute [C/»3 ’ ^ey are not 
available for computing [Cj • The simplified diagram of the computer pro¬ 
gram, Figure D. 1 of Appendix D, shows the computational order and illustrates 
the situation. This situation causes no hardship since and are not 

normally even considered in the equations of motion, In The event that terms 
such as these should be needed in some unusual case, then either numerical dif¬ 
ferentiation or iteration of the force equation will have to be resorted to. 


The gravitational force, included hefe as that due to an inverse-square 
field, can be expressed as follows: 

and J 

[ F >] * M [«aJ 


therefore. 


f< V* 

L n 


L *i 


(3.81) 


In these equations: 


[^>«] 


= the gravitational force matrix in the earth-surface 
axes, 


"M = the mass of the vehicle, 



= the sea level value of gravitational acceleration, 

= the radius of the (assumed spherical) earth, 

= the matrix representing the vector from the earth 
center to the vehicle in earth-surface axes; and 


L t» 

L .• = the column matrix of the 13, 23, and 33 elements 

*.,J ° f [*••-«]• 

The corresponding equations for the gravitational force of an oblate spheroidal 
earth are presented in Appendix C. 


Small moments due to rocket thrust can be applied to re-entry vehicles 
to provide for a lack of aerodynamic control moments at very high altitudes. 
These moments would be supplied by a reaction controller, and one form for 
such a control system is developed in Section 5. 1.7 (Equation 5.8) where the 
elements of the matrix 
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are made functions of Ot , {9 , 0 , their derivatives, and command inputs. 

The more general problem of forces and moments due to rocket thrust are 

treated in Appendix B. 
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FIGURE 3. 2 COORDINATES FOR RELATIVE MOTION OVER 
THE EARTH 
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a. Normal Orientation Angles 



b. Special Orientation Angles c. Special Orientat 


FIGURE 3. 3 ORIENTATION ANGLES 
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FIGURE 3. 4 ANGULAR DISPLACEMENTS AND CORRESPONDING 
ANGULAR VELOCITIES 
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FIGURE 3. 5 ORIENTATION OF BODY, STABILITY, AND WIND 
AXES; CONVENTIONS FOR LINEAR AND ANGULAR VELOCITIES 
AND FOR CONTROL SURFACE DEFLECTIONS 
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SECTION IV 

AERODYNAMICS AND CONTROL 
4.1 INTRODUCTION 

The basic purpose of this section is to describe methods for representing 
the aerodynamic data in the equations of motion. This description could embrace 
not only the equations for the forces and moments on the aircraft, but also the 
equations for aerodynamic heating, structural temperatures, aero-thermo-elas¬ 
ticity, control system characteristics - that is, the whole dynamic problem. How¬ 
ever, such description is considerably beyond the scope of the present work; our 
principal concern is with the forces and moments related to the rigid-body motions 
of the aircraft.. 

Before considering the equations for the forces and moments, it is neces¬ 
sary to define flight boundaries and operating limitations. Then the general aero¬ 
dynamic characteristics of re-entry vehicles can be considered as a prelude to 
formulating equations for representing the aerodynamic data. Three types of sim¬ 
ulation can be considered: (1) trajectory studies which require knowledge of just 
C L and C D ; (2) studies of rotational motions or attitude dynamics which assume 
a pre-programmed trajectory or just a constant, flight condition; (3) complete, rep¬ 
resentation for full simulation of the vehicle's flight. The significance of the ve¬ 
hicle control problem must not be missed since this represents the pilot's task - 
the real measure of what is needed in the simulation. 

Despite the voluminous publications on hypersonic vehicles and aerody¬ 
namics and despite the rapid growth of theory and experimental data, hypersonic 
flight, still remains much of a mystery and it is not possible to characterize hyper¬ 
sonic aerodynamics definitively in any general sense. Perhaps good definition will 
not be achieved until we have actually tested hypersonic aircraft, in flight. The best 
one can do is to represent the data obtained from the existing estimation methods 
even though they be inadequate, make educated guesses based on past experience, 
and hope that enough tolerance for error has been allowed. Thus, the presentation 
of numerical data has not been considered important in this study; rather, forms 
of equations have been developed which should be suitable for representing the es¬ 
timated aerodynamic characteristics of particular aircraft as they are developed. 

In this manner, the functional equations provide some measure of the complexity 
of the problem and of the requirements for simulation. 

The aircraft control problem contains many unknowns and is very much de¬ 
pendent on configuration details; at the same time the aircraft control character¬ 
istics determine to a large extent just what aerodynamic representation is needed. 
The following general criteria have been adopted (Reference 4. 5 describes detailed 
requirements): 

1. The aircraft will have an automatic control system which will 
provide the pilot with good control of the whole flight of the 
aircraft. 

2. The aircraft will have a stability augmentation system (minimum 
automatic control system) which will provide the pilot with good 
handling qualities throughout the flight: envelope. 

3. Direct manual control will be an emergency condition only, 
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and under these conditions only minimum flyable handling 
qualities can be expected. 

Pilot task analysis is the key factor - it determines the formulation of the 
problem. 

4.2 FLIGHT SPECTRUM AND CORRIDOR 

The flight spectrum can be broken into four general regions: boost, en- 
route cruise (sub-orbital, orbital or super-orbital), re-entry, and landing. Boosted 
flight may include pilot control, particularly in the event boost is aborted. The pilot's 
job during cruise depends on the task he is given, and may include the initiation of 
re-entry (the transition from cruise to re-entry). Transition is critical because 
it governs the piloting task during re-entry. Transition may be effected by firing 
retro-rockets, by aerodynamic means, or by both combined. Re-entry poses a 
piloting task under all circumstances; even with fully automatic control, the pilot 
must monitor the re-entry and be prepared to take over if a malfunction occurs or 
if unforeseen circumstances require a change in plans. The landing problem for a 
hypersonic aircraft will not be much different from that for contemporary high-speed 
aircraft, though the piloting requirements will probably be more stringent because 
speeds will be higher and less power will be available (thrust and duration for landing 
engines will be kept small because of the severe penalty incurred from excess weight). 

Re-entry, the flight region of primary concern in this study, has the purpose 
of landing the aircraft on a given spot. The range achieved following re-entry is 
very sensitive to small changes in re-entry velocity and flight path angle, so what 
would seem to be minor effects (such as earth oblateness) can take on significance. 

In order to minimize the uncertainties, it would be desirable to minimize range 
during re-entry and hence minimize dispersion. The use of maximum re-entry 
angles ( ) are thus implied, and the analysis of Lees, Hartwig and Cohen (Ref¬ 
erence 4.7) shows that angles up to - -6° are feasible. To achieve such a 

if,et will require large retro-rockets, the maximum use of drag, and possibly neg¬ 
ative lift (inverted flight). Precise control in this transition region will be required, 
including precise attitude control. Aerodynamic moments will be small, and pri¬ 
mary control will be by inertial devices or reaction jets. Though not a problem 
analyzed here in detail, it should be recognized that how the aircraft is controlled 
during this transition from cruise to re-entry will have a profound influence on re¬ 
entry. Very stringent requirements will be placed on time or spatial position, re¬ 
entry angle, attitude, and proper performance of the task. 

The "h-V diagram" is a useful and commonly used tool for representing the 
flight spectrum in re-entry (e.g. , see Figure 4. 1 in this report, also Reference 
4. 2 and Batdorf in Reference 4.6). On this daigram one can plot the various oper¬ 
ating characteristics and limit lines: equilibrium flight or glide conditions, heating 
limits, structural limits, maneuver limits, and aerodynamic limits. These define 
the corridor for safe flight. The h-V diagram is essentially a state diagram which 
defines the energy of the system - potential vs. kinetic. In addition to the h-V 
information, information is needed on the trajectory and how to control it - the 
range problem. Range control is a complex problem, especially since at high 
speeds the range is very sensitive to minor changes in atmospheric characteris¬ 
tics or aircraft characteristics and flight condition (e.g., the range is indetermi¬ 
nate if V-^Vc, and A, is high enough, i.e. , R—*-oo ). Range control is the objec¬ 
tive of what has been called "energy management"; many concepts for energy 
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management systems have been proposed and studied (e. g. , Reference 4. 3). 

ar ® not concerned with the energy management system per se, but only 
with the variables that affect it; we are not concerned with how to control range, 
ut rather with the variables needed for range control and the accurate computation 
of range given the proper control. The pertinent aspects of energy management, 
the state in h-V space and the trajectory or range in physical space are described 
very briefly in the following sections. More comprehensive treatment of energy 
management is found in References 4.3, 4.4, and 4. 13. This last paper by Yoler 
includes a concise review of the whole performance problem. 

4.2.1 Equilibrium Conditions 

Equilibrium flight assumes a balance of gravity, lift and "centrifugal force" 
the acceleration due to rotating around the earth. Equilibrium flight, a concept 
valid only for small flight path angles, is characterized by the following simple 
equations: 

L = "(t- jr) (4.i) 

where 

L = lift 

W = weight 

= vehicle geo.centric radius = earth radius 
='fcentrifugal force" 

If a non-rotating earth is assumed, then Equation 4. 1 may be written as 



and the right-hand side is completely determined for any given h and V since fi 
can be obtained from the standard atmosphere. Thus, curves of constant V//C L S 
can be plotted on the h-V diagram, or curves of constant C L can be plotted for 
a particular aircraft. Equilibrium glide curves for the study aircraft are shown 
in Figure 4. 1 for two cases: (1) C L = .7 which is characteristic of 
(2) C L for 0C = 15* which is characteristic of L /and hence the minimum 
equilibrium C L for re-entry. 

If earth rotation is included, then Equation 4.2 becomes 

= A#yL ( 4 . 3 ) 

C L S I - V 

9 r 

and may be computed to sufficient accuracy from the relation 

V* * vj + Ver + 2Ver V* sin r (4.4) 

where 

= azimuth angle of flight path 

Vir = velocity due to earth rotation at the 
particular latitude, -O-r ^ 


58 



NAVTRADEVCEN 801A 
Equation 4.4 may be solved for y r to yield 

14 - - Kir t T *\f V*- »4r* (4. 5) 


and the + sign is used except when both ^ 180* and V A < - V-., SHiflf. Since 
Ver is a function of latitude, it is clear that a single h-V diagram cannot account 
properly for the difference between inertial and relative velocity except in the case 
of a nearly equatorial re-entry. However, by using an average value of I4r » an 
acceptably accurate picture will usually be obtained. An upper bound for the error 
is shown in Figure 4. 1 by equilibrium glide curve for an eastward equatorial 
re-entry. 

The range for an equilibrium glide can be roughly approximated by the 
equation 

* r &)' U (l-$ r ) <4 ' 6> 

where V Q = initial velocity at start of glide. 

It is seen that while C L controls the profile in the h-V plane, it is L/D that con¬ 
trols the range. Thus, a high L/D ma x is needed in order to have a large range, 
and control of L/D is needed in order to control the range or to vary the landing 
spot. 


4.2.2 Aerodynamic Heating 

To a great extent the aerodynamic heating depends on the aircraft shape 
and physical characteristics: the shape and thickness of the nose, wing and tail 
leading edges, etc.; the range of 0* and /f allowed in flight; the structure and 
method of cooling it. The basic problem in flight is to keep the temperatures at 
various critical points on the aircraft less than the allowable limit temperatures. 
Various temperature-limit lines ( T umr lines) can be plotted on the h-V dia¬ 
gram, and these will usually depend on 0* or C L . The calculation of T UM/r 
lines is complicated; however, theoretical heat transfer rates for simple shapes 
and conditions do give an insight to the problem. Lees (Reference 4. 12) gives an 
equation developed from the Fay-Riddell results for the stagnation-point heat- 
transfer rate (convective): 

<5, (4.7) 


where 


C - constant depending on the characteristics of 
the air and the shape of the surface 

R = radius of curvature of the surface. 


A curve for Q f = 200 Btu/sec-ft 2 has been plotted on the h-V diagram in Figure 
4. 1 to typify a heating limit curve for the study vehicle. This curve represents 
for the nose of the vehicle which has been assumed to have a radius 
R - .5 feet 


and the constant, C , according to Lees, is 
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/* - is 

C - * n EARTH 

= 1 for axially symmetric flow 
Dearth - 15.5. 


Of course, skin and structural temperatures depend on more than just heat trans¬ 
fer rates. Not only must convective, conductive and radiative heat transfer rates 
be known, but also the heat flux through the structure must be known, and the ef¬ 
fect of cooling by ablation or other methods must be included in order to solve for 
the actual temperatures. 


4.2.3 Structural Limits 

A number of limitations fall in the structural category: air loads, maneuver 
loads, flutter (particularly panel flutter). The structural problem has been sum¬ 
marized by Batdorf in Reference 4. 6, and References 4. 1, 4.2 and 4.4 contain 
analyses of the maneuver loads problem. The primary influence on the structural 
limits is the dynamic pressure ), so that lines of constant (f.[ in the h-V plane 
are typical of structural limit lines. A curve for = 500 psf has been plotted 
on Figure 4. 1 and is typical of a structural limit line for a re-entry glider with 
low wing loading. For an actual aircraft the structural limit lines, as in the case 
of T umit lines, depend on the vehicle's attitude and flight condition. In particular, 
heating is significant since high temperatures and high stress combine to give crit¬ 
ical conditions. 


4.2.4 Maneuver Limits 

The maneuver limits for a re-entry vehicle are similar to those of a con¬ 
ventional aircraft (load factor, or "g" limits), but an additional curve appears on 
the h-V diagram called the "recovery ceiling". A recovery ceiling, as it might 
appear, is sketched on the h-V diagram of Figure 4. 1. The recovery ceiling is 
defined as the maximum allowable altitude (potential energy) for horizontal flight 
(apogee) at a particular velocity. If this altitude is exceeded, then on re-entry 
the vehicle will exceed one of the flight boundaries - burn up, break up, or crash 
into the ground. This will be true if the recovery ceiling is exceeded for the given 
velocity, irrespective of the direction of flight (i.e. , horizontal flight is only the 
least critical case). The recovery ceiling is primarily of interest to the boost phase 
of hypersonic flight. It is desirable to avoid crossing this boundary at any time 
during boost, so that abortion and subsequent re-entry are always safe. 

Because of the low values of L/D in hypersonic flight, drag is an impor¬ 
tant factor in the maneuver loads. The maximum load factor (maximum total "g") 
on the h-V diagram is essentially a line of constant for any given angle of 
attack. Maneuver load boundaries are discussed further in References 4.2 and 
4.4. 


4. 2. 5 Stability and Control 

The over-all problem of stability and control is a complex one, involving 
handling qualities under various types of operation (automatic control, manual 
control with stability augmentation, direct manual control), maneuver limits 
(accelerations on the pilot), attitude and control limits (stability, heating, and 
control authority), and range control (energy management). The stability and con¬ 
trol problem is intimately linked with the whole flight spectrum, and what 
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limitations will be imposed by these problems are not presently clear. This is, 
in fact, the purpose behind the present study: to develop equations which can be 
used to help define the critical stability and control problems. However, certain 
characteristics are fundamental and can be described with present knowledge. 


The basic problem is that of providing control of angle of attack ( ) and 

bank angle ( 0 ) in such a way as to satisfy the handling qualities requirements, 
performance requirements, and heating limitations, and do this for each of the 
control modes (automatic, manual, etc.). The angle of attack range during re¬ 
entry is limited at the low end by 0C for L/D max > at the high end b Y stability and 
control limitations; contrary to normal flight, re-entry flight occurs essentially 
on the "back-side of the power-required curve". The reason is that for 
heating increases drastically (Reference 4. 3). It is desirable to have 0C range 
from to ’ and from Figure 4.3, this corresponds to roughly 

15* <■ < 50*. It would be of value to obtain maximum C D with L/D = 0 (i.e. , 

= 90°), but it is unlikely that the severe stability and control problems pre¬ 
sented by flight at 0t = 90* would justify it. Achieving & = 50* is difficult enough. 

Range and altitude control is achieved principally through bank angle ( 0 ). 
Though the equilibrium glide is continually discussed as the mode for flight during 
re-entry, it is a condition which is virtually impossible to achieve except for a 
non-lifting vehicle. The actual flight path achieved involves skipping (See Figure 
5.3 for example). The first skip is essentially impossible to avoid (the first over¬ 
shoot below equilibrium glide), but the best method of damping the skip is by 
banking. That is, 0 is used to reduce "lift" (normal to the earth) while main¬ 
taining drag. Thus 0 becomes effectively a control for U and L?/D independent 
of & ( L? denotes the effective lift normal to the earth). In addition to its function 
as an altitude control, 0 controls the lateral range. 

The manner in which & and 0 control the over-all aircraft motions can 
be seen from the following simplified equations of motion in which the effects of 
earth rotation (among others) have been neglected. 

L cos 0 -- W(/- -H.V* (4-8) 

L s in 0 = ->* Vf T (4.9) 

>nV - 0 (4.10) 


where 


if - the flight path elevation angle, assumed small 


The underlined term 
The range, from the 


= the flight heading angle. 

) is zero for equilibrium glide by definition ( 
initial conditions, is given approximately by 


x */.' v«. r T at 


Y - f** V S'n ^7. d.t 


y = o). 

(4.11) 


where X and Y are the longitudinal and lateral range, respectively. Now it is 
clear that L cos 0 takes the place of L in the equilibrium glide equation (Equa¬ 
tion 4. 1). L/D max still produces the maximum range, but the total range is re¬ 
duced by banking. It can be shown that 0 = 45* produces approximately the 
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maximum lateral range. By banking successively left and right while holding a 
high 0C , the longitudinal range can be decreased without changing heading. The 
effect of 0 control on range is sketched in Figure 4.2 (including two digital so¬ 
lutions for the study vehicle), and the envelope of possible landing points (called 
a foot print) is indicated for the particular initial re-entry condition. 

The method of controlling the skip - of damping the hypersonic phugoid - 
is readily apparent. L cos 0 controls y (Equation 4.8), and so once if = 0 
at the bottom of the first skip then 0 can be controlled to hold if 0 and main¬ 
tain if =0. Thus, Of can be used to control range, through D, and 0 used to 
control h through L cos 0 , Such control takes on particular significance for 
super-orbital re-entry ( V ^ V c ). Here negative "lift" ( L* ) is needed to stay 
in the atmosphere (Reference 4. 7) : and inverted flight ( 0 = 180*) is the only 
reasonable method of achieving it since the use of negative Ot creates many prob¬ 
lems (heating, structural, negative "g" on the pilot). For sub-orbital re-entries 
0 < 90* would normally prevail; though if the aircraft were permitted to start 
skipping out, then inverted flight ( 0 - 180*) would be required. Such a tech¬ 
nique might be used to gain altitude. Roll rates ( 0 ) for controlling range might 
be quite low, but high rates (as measured by the required aerodynamic or aileron 
control power) may be needed for controlling the skip. 

4.2.6 Pilot's Displays 

Certain quantities ( A. , V , 0C , etc.) take on particular significance in 
terms of the piloting task in re-entry. Once initiated, the piloting task in re¬ 
entry is essentially that of controlling the trajectory through the attitude of the 
aircraft while staying within a limited h-V corridor. The pilot must be supplied 
with information on range, the aircraft's attitude, and the corridor (output, input, 
and constraints). The essential features of this control problem, the energy man¬ 
agement problem, have already been summarized briefly using simplified equa¬ 
tions and concepts to demonstrate the gross effects. There is considerable liter¬ 
ature on energy management; Reference 4.3 (Unclassified) gives a good compre¬ 
hensive treatment of the problem and includes many detailed results. Three pri¬ 
mary displays are suggested in Reference 4. 3: a range display showing the es¬ 
timated footprint and the vehicle track on it, an h-V display with the various limits 
and the predicted path in h-V space, and an attitude display similar to the one pre¬ 
sently used but with sensitivities increased where precise attitude control is re¬ 
quired. The problem of how to measure, process, and compute the necessary 
quantities should not be minimized; it is an extremely difficult one. However, 
the problem for this study is a different one; it is the problem of computing from 
the equations of motion the quantities necessary to simulate the displays. It is 
logical to include Oi and in the primary display, perhaps also y and T . In 
addition to the normal aircraft instruments, there will be a need for instruments 
associated with special problems and equipment (i. e. , energy management system, 
propulsion system, reaction controls, life support systems, etc.). Thus a sim¬ 
ulator will have to include not only the equations of motion for the vehicle, but 
also equations for computing all the displayed quantities for the subsystems and 
any special sensors required (e.g., ones to measure skin temperatures at criti¬ 
cal points on the vehicle). The equations of motion developed in this report have 


* An interesting possibility is suggested. By flying at escape velocity (36,000 
fps) in inverted flight, a normal one-g environment would be achieved for per-, 
sonnel in the airplane. 
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been derived with these facts in mind. Especially, the axis systems and vari¬ 
ables used or solved for in the equations were selected to provide the required 
data as expeditiously as possible. 

4.3 GENERAL AERODYNAMICS FOR RE-ENTRY VEHICLE 


The aerodynamic problem for a re-entry vehicle is essentially that of de¬ 
termining the aerodynamic forces and moments acting on the vehicle as a func¬ 
tion of the parameters (flight condition, attitude, and shape) that affect these 
forces and moments. It is not possible to give explicit functions for the force 
and moments except for a particular vehicle. However, functional forms can be 
given, and to this intent we start with some definitions and assumptions. We as¬ 
sign an axis system ( X , y , % ) with specified origin, and then write the forces 
and moments in terms of the air density, velocity, and the vehicle's wing area, 
span, and chord_. 


M * t /> v i * [<v] 




[c,]= C, , [e„]a ci , M= [„ c d 


Isj (4. 12) 


F r and M t are reaction forces and moments due to thrust and reaction 
controls which are not conveniently expressed in coefficient form. 

The forces and moments can be sub-divided into certain groups with which cer¬ 
tain parameters can be associated. Without loss of generality, we may write 
the forces and moments in differential form: 
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The notation in Equation 4. 13 is much abbreviated and implies the following 
concept. The coefficients of the variables are in general functions of the var¬ 
iable in question and all preceding variables £e.g. , s with all 

remaining variables = 0, but *e» ^ « *S r t -fik, > , 

The terms and S <*i. refer to structural deformations in the i th mode due 
to loads and the i^ mode due to heating, respectively. A similar expression 
can be written for eLC-^ . With the defined meaning of the notation in Equation 
4. 13, the forces and moments at any time can be obtained by performing the 
multiple integration 

r. i 

(4. 14) 


[ d( *] * N 


where the constant of integration is a function of Mach number alone 

K*. 


N- 


C Vt CM) 
Cm CM) 




*• " J ...«o 

If , C- , are zero (as is usually assumed), then the \C-f\ integration 

(Equation 4. 14) can be performed independently of time (time dependence expres¬ 
sible in terms of dependent variables and their time derivatives) and yields a 
general expression for the force and moment coefficients. The constant [&/$] 
is evaluated starting from M = 0 and some convenient Reynolds number ( R m , ). 
The above form of the equations allows for non-linear variations of any of the 
derivatives with respect to M , , dt , /S , . .., but evaluation would be a fear¬ 

some thing. In practice, drastic simplifications are rrade. 


With the present knowledge and methods of estimating aerodynamic de¬ 
rivatives, and based on experience, we can approximate the aerodynamic forces 
and moments in the following manner. There are basically four regions of flight 
defined by the flight velocity : 


Subsonic M < 1 

Transonic MSI 

Supersonic M > 1 

Hypersonic M ^ 1 


In the first three, the over-all characteristics tend to be linear and to depend on 
Mach number. For hypersonic flight, the over-all characteristics tend to be 
non-linear and to depend on Reynolds number and aerodynamic heating. Accord- 
ingly, we find it convenient to treat this problem in two parts. First, we con¬ 
struct the characteristics for 

M <■ Ml 

where Ml i8 the Mach number at which hypersonic effects become dominant for 
the particular vehicle involved (typically 5 < Ml 15). Next, we add the effects 
for M > Ml. 


In the lower speed regions, M < Ml , we make the following assumptions, 
implying a well-behaved airplane: 

1. <t will be small, less than stall 

2. fi will be small, less than a few degrees. 

Under these assumptions, the forces and moments can be expressed as follows. 
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Equations 4. 16 and 4. 17 include the terms that will normally be encountered; 
additional terms which may be needed in special cases are discussed below. 

The notation implies incrementals (in theory resulting from each successive 
integration in Equation 4. 14). For non-linear terms, the increment is attri¬ 
buted to the first variable that it is said to be a function of [e.g.. (*,M) 

is the increment in C m due to Ot with 0 , S e , S* , S r , etc. all zero, and it 
is a function of M and non-linear in (Kl. For linear terms the increment is 
attributed to the variable with respect to which the partial differentiation is 
taken. Aeroelasticity may be accounted for in the derivatives (quasi-static ef¬ 
fects). Control interaction between ailerons and elevator (common for delta 
wings) is included. If configuration changes (e.g. , flaps and gear) have signify 
icant effects, appropriate terms must be added. For large sideslip angles or 
aileron deflections, pitching moment increments may be required [i. e. , 

(/*»*,*) and AC m (S a , at, M), respectively J , and Cj and C* may 
show non-linearity with/tf or . Except for unusual effects peculiar to a par¬ 
ticular airplane (e.g. , wing dropping at transonic speeds). Equations 4. 16 and 
4. 17 should be adequate to cope with all significant aerodynamic data for M 


In the hypersonic region, the non-linear effects associated with 

hypersonic flow become important (References 4.8 and 4.9 describe the gen¬ 
eral characteristics). The elements of the phenomena may be shown as follows. 
The linearized theory for supersonic flow((given for example in Reference 4. 10) 
gives for the pressure coefficient 




(4. 18) 


where. 0 is the flow deflection angle from free stream (positive if the normal 
to the surface points upstream). The simple Newtonian expression for the hyper¬ 
sonic pressure coefficient is given as 


Cf, « Z tin* a 


(4. 19) 


f 
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where $ is the angle the surface makes with the oncoming flow (Newtonian 
theory does not specify flow on the back surface - ambient or less). Consider 
what happens as Mach number becomes large: 

BECOM0S VAL,t> 

The same phenomena is demonstrated for the lift on a two-dimensional flat plate: . 


C L (supersonic) = ( 4 . 20 ) 

! M-l 

C L (HYPERSONIC) - CM a, • 2*i* Z 0c 6tt & (4.21) 


One can see that at high Mach numbers Equation 4.21 effectively adds a correc¬ 
tion to Equation 4. 20; that is, the sum of the two terms would give the correct ans¬ 
wer. Furthermore, for small CC the sum of Equations 4.20 and 4.21 would also 
give the correct answer at supersonic speeds. This formulation suggests that ap¬ 
propriate expressions for the hypersonic forces and moments can be obtained by 
adding corrective terms to the expressions for the supersonic forces and moments. 
This, in fact, does prove to be a convenient method of handling the data; but the 
problem is considerably more complicated than Equations 4. 18 through 4.21 indicate. 

The significant factors affecting the representation of hypersonic aerody¬ 
namics may be summarized briefly as follows. 

Mach Number Independence: The changes in the aerodynamic coefficients 
with Mach number are very small and can generally be neglected. 

Viscous Effects : Reynolds numbers are low, boundary layers are large, 
and there is strong shock-wave boundary-layer interaction. As a result the aero¬ 
dynamic coefficients and particularly the drag coefficient are dependent on altitude. 
The effect is demonstrated by the data in Figure 4.4. 

Longitudinal Trim: The data of Figure 4. 3 shows that QL —►SO* for C L 
At hypersonic speeds the center of pressure is near the centroid of area, mud?** 
farther aft than at subsonic speeds so that large nose-up moments will be required 
for trim. These may be difficult to achieve without variable geometry. The most 
efficient method would be to extend a nose, surface; a simpler w^y would be to tip 
up the nose. Conversely, the basic configuration could be designed to trim at high 
angles of attack at hypersonic speeds, and then tail surfaces could be extended 
(perhaps folding out from vertical fins) at lower speeds to provide the needed sta¬ 
bility. The primary effect of the geometric change would be on C* , Cj. , and C ^ 
but other derivatives such as might also be strongly affected. ' 

Non-linearities in Statics: As evidenced by the data of Figure 4.3, it will 
be necessary to represent all the forces and moments as function of angle of at¬ 
tack. However, it is clear that any lifting surface will have strongly non-linear 
properties. Accordingly, the forces and moments will all have to be represented 
as functions of the two attitude variables 0C and over the expected range of 
these variables. 

Non-linearities in Control Effectiveness: The forces and moments due to 
control surface deflection will be found strongly non-linear. For example, in the 
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study configuration (see Figure 5. 1) the elevator effectiveness is small for 

>■ & (shielded by the wing), but increases rapidly as the elevator is de- 
fleeted down, - S 6 <■ & . The rolling and yawing moments due to aileron deflec¬ 
tion (differential elevator deflection, ) will be non-linear with 8 ^ and will 
depend strongly on 0C and S e and less strongly on /S The vertical tails will 
be ineffective at small angles of sideslip unless the rudders are deflected differ¬ 
entially outboard ( fr L > 0; % r . m <. 0) in their neutral position, and it may be ex¬ 
pected that only outboard rudder motions will be used (for a single tail, a split 
trailing edge, opening to a wedge at hypersonic speeds, could give the same ef¬ 
fect), Thus the rolling and yawing moments due to rudder deflection ( 8 r ) will 
be strongly non-linear with 8 r , 0t , and /B . What is true for the study ve¬ 
hicle is typical of what can be expected: control effectiveness will be non-linear, 
will be affected by attitude ( & and fi ), and there will be strong interactions of 
the controls; furthermore the static stability of the vehicle ( C >*. , , C* ) 

will also be a function of control surface deflections. * * 


Rotary Derivatives : Only small values of the non-dimensional or reduced 
rotation rates {job/zv , £e/,2|/ , ... &b/Z V ) will be encountered at hyper¬ 
sonic speeds, and aerodynamic damping will be minimal. Etkin has shown this to 
be true for the longitudinal short period mode (Reference 1. 6), and the same can 
be assumed for lateral or coupled motions. Accordingly, only rough approxima¬ 
tions for the rotary derivatives ( , • • • , etc.) will be needed. 

Certainly they can be linearized with respect to rotation and at most made functions 
of Ot and /3 (just OC if is not too large). 




Heating Effects: Substantial structural deflections may be caused by aero- 
dynamic heating, and if so, then the effect of such deformations will have to be 
included in the aerodynamic forces and moments. It can be surmised that only the 
moments due to such deformations will be significant. It may not be sufficient to 
consider equilibrium or steady state conditions, more likely, it will be necessary 
to include degrees of freedom associated with heating. The significance of defor¬ 
mations due to heating will depend on the particular configuration; but, in general, 
if such deformations are significant it will be necessary to compute transient tem¬ 
peratures, then the deformations, and finally the corresponding aerodynamic 
moments. 


From the foregoing discussion of the various factors affecting the repre¬ 
sentation of hypersonic aerodynamics, it is clear that most of the problems are 
intimately associated with the specific vehicle configuration employed. However, 
one can formulate the likely functional dependencies by making certain assump¬ 
tions, and thus indicate the problems which are likely to be encountered. With 
the study configuration in mind (shown in Figure 5. 1), the following assumptions 
are made to arrive at an aerodynamic representation for hypersonic flight. 

1. Mach number independence 

2. Viscous and dissociation effects depend on altitude alone 

3. Variable nose tip-up is used to trim at high Ot 

4. Rudders are deflected outboard (towed out) to maintain 
directional stability and rudder effectiveness 

5. Moderate sideslip angles ( ) 

6. Heating is symmetric and produces only symmetric deflections 
of the structure (good only for small ^ ) and these deflections 
can be characterized by modes. 


2 
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The following notation is adopted: 

nose tip-up deflection 

differential antisymmetric rudder deflection 
deflection of the i*^ mode due to heating 

The aerodynamic force and moment coefficients for hypersonic speeds can now 
be written in functional form, where terms and functional dependencies that 
are likely to be of lesser significance are omitted. 
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(4.23) 

Again, the notation ia such that the first variable in the incremental term is the 
one to which the effect is attributed [ i. e., 4 Cj (S^, 0, 0) mean the change in 
C 4 due to , and this is non-linear with % 0 and is a function of 0 and 

«]. It should be noted that the conventional derivatives may have contribu¬ 
tions from many sources. For example, it is implied that • 9 will 

be affected by 0 , Ot , S* , 5 # . S r , , and % m0 . 'Hie representation 

will require, in some cases, generating functions of four variables. The 
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4s,. must come from equations determined by a detailed analysis of the aero-thermo¬ 
elastic properties of the vehicle. 

Having developed one general formulation for the aerodynamic forces at 
subsonic, transonic, and supersonic speeds and a second general formulation for 
hypersonic speeds, we may well ask how does one combine the two? The super¬ 
sonic derivatives, which are needed only for moderate angle of attack, should be 
asymptotic to the hypersonic values at some low altitude, say between 150,000 and 
200,000 feet where the Reynolds number is still large. Thus, the supersonic values 
can be used as the basis for hypersonic forces and moments, and the additional ef¬ 
fects associated with hypersonic flight can be added to the supersonic formulation 
for (it is presumed that at altitudes above 200,0Q0 feet where h 

begins to have an effect). Where such a scheme is not possible because of pecul¬ 
iarities in the flight profile, it may be necessary to make the formulation continuous 
with Mach number for all M. This is discussed in further detail in Section 4.4. 


4. 4 METHODS OF HANDLING AERODYNAMIC DATA 


In this section the actual methods for handling the aerodynamic data rela¬ 
tive to the functional dependencies are discussed, and equations and numerical 
data are presented from which some appreciation of the magnitudes of the various 
effects can be obtained. Then the data required for various types of simulation are 
discussed. In particular, three types of simulation are considered in which the 
requirements are formulated by the task being simulated: trajectory control, at¬ 
titude control, complete simulation. 


c 


4.4. 1 Aerodynamic Data 

As pointed out in Section 4.3, the aerodynamic data for flight at subsonic, 
transonic, and low supersonic speeds is essentially linear and depends principally 
on M and Qt . As velocity increases to hypersonic speeds, non-linearity pre¬ 
dominates and the stability derivative concept is no longer valid. However, for 
very small disturbances (in ft , ^ , S # , etc.) the aerodynamic characteristics 

at hypersonic speeds are similar to those for lower speeds, and the values of the 
supersonic stability derivatives approach the hypersonic values in somewhat of an 
asymptotic fashion. The Vetter scale (after Hans Vetter of Douglas Aircraft who 
originated it) provides a very useful way of plotting the data so that the pertinent 
aspects of the various flow regimes are put in proper perspective. The Vetter 
scale is given by the following formula: 

1= M 0<- /*«*/ 

Z • 2 - pf t<./1*.oo 

where X is the abscissa. The Vetter scale compresses the whole Mach number 
range into a finite scale, and its utility becomes apparent when one considers the 
following characteristics. 

Prandtl-Glauert Rule: 

,4 ' 24 ’ 

Linearized Two-Dimensional Supersonic Wing Theory: 


<M.,r 


vb m £['+*(*?]*• 


(4. 25) 
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Data for the normal force as a function of M , Ot , and A are shown in 
Figure 4.4 for a typical hypersonic glider and serve to illustrate most of the 
problems that are encountered. Part a shows that the Vetter plot of is a 

nicely behaved function. Part b shows the effect of Ot and . The presence 

of hypersonic effects is dependent on Ot ; for large Ot , hypersonic type flow 
occurs at lower H than for small Ot (the familiar effect of shape on the exis¬ 
tence of hypersonic flow). In fact, the behavior of the forces on a wing in stalled 
flight at subsonic speeds is somewhat like the behavior at hypersonic speeds. 

Part c: shows the interaction of viscous effects and shock effects (including dis¬ 
sociation) which causes the forces to be altitude-dependent. Though the data shown 
in Figure 4.4 characterizes the behavior of Cj.* in symmetrical flight, it also 
characterizes the behavior of most of the other small disturbance derivatives 
( e -g. i C* , C n ^ , ) since they are either closely related to C (e.g. , 

is*proportional to^ ) or are an analogous derivative (e.g. , c* is 
the analogous lateral derivative). ’ * 

To further characterize the behavior of the hypersonic forces and moments, 
we examine the normal force on a flat plate and a cone. For the flat plate, the 
normal force is given by 

<hi ~ *($) (4.26) 

where is the normal component of velocity. In terms of C* as a function of 
Ot and /S , Equation 4.26 becomes 

~ c i 2 (f)\r\ 

For a cone with half-cone angle Sy 
ence area) is given for Ot < S y as 

C N * cos* Sy sin 2 at (4.28) 

Then Cy and Cj are found from Equation 4. 28 to be 

-C 7 = * *y (f)(f) x 2 c*‘ *y tin ot os 0L~* V (4. 29) 


Cy = 2 cct 2 2 cot* $ v tin/S cot ^3 c»t Ot (4.30) 


For small disturbances about Ot - /S = O , the force derivatives are as follows: 

Flat Plate: 


= 0 

_ *0L 

Cone: 

(4.31) 

“ C £ = 2 cos 2 $ y 

(4.32) 

— Z COS 2 5y 

(4. 33) 

Equations 4.26 through 4.33 demonstrate 
Whereas at Ot - 0 a flat plate will have C L = 0, 

several interesting properties, 
a slender body will have a 


« Z si» ot ctfi | si* * coo/3 | (4. 27) 

, the normal force (with base area as refer- 
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Unite C(. A Furthermore, it is shown in Reference 4.8 {Chapter 3, page 82) 
that Equations 4. 32 and 4. 33 agree with the exact gas dynamic solutions for 
cone lift at small angles of attack. These results hold for 5*< S v < 30* and 
2 < 20 (presumably also for M ^ 20). Thus, the supersonic small dis¬ 

turbance derivatives are shown to be asymptotic to the hypersonic values for 
both a cone and a flat plate. The gross behavior of wings and bodies for large 
Ot and ^ described by Equations 4.26 through 4.30; it is clear that the direc¬ 
tion cosines ic/v, sv~/V , and c*r/V are more natural coordinates than 01 and /6 
for describing the forces and moments. 

It is instructive to look at flap effectiveness from the Newtonian view¬ 
point. For a flat plate with a flap, the normal force on the flap is given as 

Cfl * 2 *;»*(& + S), As s 
C„' = o ,A>* S >-* 

Equation 4. 34 expresses the effect of wing shielding of the flap. The total nor¬ 
mal force on the flat plate, where Sf is the flap area, is given by 

Cry — 2 s in 2 Oc j ^S/Vt * t ^ 

Equation 4.35 shows that both the lift curve slope and flap effectiveness are af¬ 
fected by flap deflection. The results can be extended readily (in the same manner 
as for previous equations) to determine control surface derivatives, and partic¬ 
ularly to show the necessity for antisymmetric rudder deflection ( ) if a sta¬ 
ble is to be obtained. 

4.4.2 Various Applications 

In a flight simulator the objectives govern the degree of sophistication re¬ 
quired in the aerodynamic representation. If only energy management is of inter¬ 
est - not the attitude control problem - then adequate simulation may be achieved 
by just representing C L and C 0 (or C * and C» ) as a function of dt , At and 
'H.. Some nominal transfer function could be used between the pilot's control stick 
and Oi . If C L and are to correspond correctly with <L , then the effect of trim 
should be included. 

If only the attitude control problem is being simulated, two possible sim¬ 
ulation techniques are suggested: one with constant speed and altitude (fixed flight 
condition), the other with a time programmed speed and altitude. In either event 
the forces and moments would have to be computed as functions of the attitude and 
control deflection ( 0C , £ , ^ , r , , etc.); but, in the first 

case no change in flight condition would be called for, in the second case the only 
additional variable would be t (i.e., the programmed M(t) , Aft) , etc. 

could be used to eliminate any explicit dependence of the forces and moments on 
flight condition). 

If complete simulation of both the attitude and trajectory control problem 
is required, then the full aerodynamic representation will be required. However, 
how complicated this has to be depends to a great extent on the leeway allowed 
the pilot. If he is required to follow a particular schedule of V , £ , and 0C quite 


(4. 34) 
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closely and hold lateral perturbations to a minimum {fi , -f> . and r small), 
then a linear representation (about a trim a ) using just stability derivatives 
might be adequate. On the other hand, if large perturbations are required in the 
simulation (the ultimate might include no restrictions on attitude at all, so that 
tumbling or spinning could be simulated), then the non-linearities associated 
with these perturbations must be represented. It may be presumed that if the 
pilot is flying with direct manual control (probably only an emergency condition), 
then large motions will occur - mostly at the higher altitudes. In this case, ac¬ 
curate simulation will require the full non-linear representation of the aerodynamics. 

4. 5 VEHICLE CONTROL 

Though the control problem in re-entry and landing involves both trajec¬ 
tory and attitude control, our concern here is only with the latter. It has been 
shown in Section 4. 2 that the pilot's task is essentially to control a. and 0 , both 
through wide ranges. We may add /S to the list; but the pilot will want /S — 0. 

This is essentially the same task that the pilot of a conventional aircraft performs, 
but with different emphasis. For one, much larger 's will be used in re-entry; 
two, the in re-entry will be much smaller. In fact, normal acceleration - 

which is so important an input to the pilot in conventional flight - will be almost 
non-existent ( S' Q ) in the initial part of the re-entry. 

The details of the actual flight control systems that will be used in re-entry 
aircraft are not of real concern here - they depend too much on the specific con¬ 
figuration and hardware involved. But the aircraft response that the flight, control 
system must produce , if the aircraft is to have satisfactory handling qualities, is 
not subject to radical variation (see Reference 4. 5). Since our concern is with the 
piloting task, we need not consider regulatory control modes (i.e. , attitude hold 
mode, altitude hold mode, etc.) but only actual responses to pilot inputs. 

Two types of control systems can be envisioned: direct manual control and 
control through an automatic stabilization system. Two types of control actuators 
are involved: aerodynamic control surfaces and reaction controls. Direct manual 
control needs no further discussion, except perhaps to state that its implementation 
is straightforward, and to restate that it will probably present a severe control 
task for the pilot. With automatic stabilization, it will make little difference if 
aerodynamic or reaction controls are used if a desirable response is produced. 
Differences will arise from the failure to produce such a response, perhaps due to 
control limits, limit cycles, unforeseen problems, etc. Thus at this stage the 
question becomes what is a desirable response. 

From the handling qualities research (summarized in Reference 4.5) there 
seems to be little doubt as to the formulation of a desirable response in angle of 
attack. Pilots like the following transfer function. 



f 2 0.7 co n =0.7 cps 

The value of K is not well established and admits to considerable variation. There 
is every reason to expect the above criteria to hold for re-entry. 

Much less is known about lateral handling qualities - perhaps because they 
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are more complicated than longitudinal ones. However, the response of and 
Y* should be well damped. Two possibilities for 0 are suggested: 



JL = K 

(4.37) 

and 

_0_- * 

*S, s* 21 

(4.38) 


Equation 4.37 is the conventional response expected for a lateral stick (aileron) 
deflection; but for re-entry Equation 4. 38 has much to recommend it. Perhaps 
some combination will prove best. 


In the computer studies performed for this report, a perfectly adaptive auto¬ 
matic control system was assumed - one that would produce the desired response 
irrespective of the aerodynamic moments. Such a control system would be capable 
of producing the following moment equations. 

L » 10 - 0 c)* & 

* 

v * o -a) < 4 - 39 > 

Then the transfer functions for command inputs and small disturbances from 
straight symmetric flight would be the following 

cc _ _ 

^_ 

A “>>> s * ***1 

0 __ ^0 _ 

4 (4.40) 


where 

c*>rv a s ^ Ay 2 

"v = * v 

= Ai,/ r x . Ar 


(4.41) 


and 

^ Ay * ^ Z-m. * * ” "A 

^4 = *** * *V 

^ A* ~ A* (4. 42) 
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Equation 4 39 was used in the computer studies and provided adequate control, 
held fivTI T 8 ? ad * to a da Pt**\t to produce a constant & ; all gains were 

t In * rtia couphng effects were encountered in rolling maneuvers (see 
i 8 Very sma11 steady state errors persisted in dr , A , and 0 due to 
a P„J.!rfl”5 t 0n8 t nd an 8 ular velocities associated with earth rotation and flying 
was imm flight path around the earth. Also, a surprisingly large transient effect 
^u red durin 8 rolling maneuvers which was associated with earth-rotation 
terms in the equations for a and ^ (see Figure 6. 1 and discussion in Section 
* ' ortunateiy, there was insufficient time to make a thorough study of the 

implications of the effect. 


4.6 REFERENCES 

4. 1 Eggers, Alfred J. , Jr. : 'The Possibility of a Safe Landing." Chapter 13. 
Space Technology, Howard S. Seifert Editor, John Wiley and Sons, Inc. , 
New York, 1959. 


Metcalf, H.: Some Consideration of Shape and Control for Hypersonic 
Vehicles. IAS Paper No. 59-144, October 1959. 

General Electric Company: Flight Control Study of a Manned Re-entry 
Vehicle. WADD Technical Report 60-695, July I960. 

Eggleston, John M. and John W. Young: Trajectory Control for Vehicles 
Entering the Earth's Atmosphere at Small Flight-Path Angles. NASA 
Memorandum 1-19-59L, February 1959. 

Breuhaus, W.O., P.A. Reynolds, and E. A.Kidd: Handling Qualities 
Requirements for Hyper-Velocity Aircraft. (Title Unclassified). Cor¬ 
nell Aeronautical Laboratory, Inc., Report No. TC-1332-F-1, September 
30, 1959, Revised January 28, I960. CONFIDENTIAL. 


4.6 Eggers: "Performance of Long Range Hypervelocity Vehicles". 

Hansen: "Some Characteristics of the Upper Atmosphere Pertaining 
to Hypervelocity Flight". 

Batdorf: "Structural Problems in Hypersonic Flight". 

Lees: "Recent Developments in Hypersonic Flow". 

Stalder: "A Survey of Heat Transfer Problems Encountered by 
Hypersonic Aircraft". 


74 



NAVTRADEVCEN 801A 


Reference 4.6 (continued) 

Ferri, Feldman, Daskin: "The Use of Lift for Re-Entry from 
Satellite Trajectories". 

Hypervelocity Flight Issue, Jet Propulsion, November 1957. 

4.7 Lees, Lester, Frederic W. Hartwig and Clarence B. Cohen: "Use of 
Aerodynamic Lift During Entry into Earth's Atmosphere". ARS Journal, 
September 1959. 

4.8 Truitt, Robert W. : Hypersonic Aerodynamics. The Ronald Press Co. , 
New York, 1959. 

4.9 Hayes, Wallace D. and Ronald F. Probstein: Hypersonic Flow Theory. 
Academic Press, New York, 1959. 

4.10 Kuethe, A.M. andJ.D. Schetzer: Foundations of Aerodynamics. 

John Wiley and Sons, Inc., New York, 1950 Chapter 11. 

4.11 Hoerner, Sighard F. : Fluid-Dynamic Drag. Published by Author, 1958 
Chapter 18. 

4. 12 Lees, Lester: "Recovery Dynamics - Heat Transfer at Hypersonic Speeds 
in a Planetary Atmosphere". Chapter 12, Space Technology, Howards. 
Seifert Editor, John Wiley and Sons, Inc. , New York, 1959. 

4. 13 Yoler, Yusuf A. : "Dyna-Soar: A Review of the Technology." Aerospace 
Engineering, August 1961. 


75 



NAV TRADE VC EN 801A 


ALL CURVES COMPUTED ASSUMING NONROTATING 
EARTH EXCEPT FOR ONE AS MARKED 


DENSITY FROM STANDARD ATMOSPHERE 



FIGURE 4. 1 h - V DIAGRAM FOR THE STUDY VEHICLE 
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INITIAL CONDITIONS: Projection Technique: Solutions 

V = 24,000 ft/ sec plotted on rectangular A - ^ 

>C = 300,000 ft grid - range computed using 

# = -1 • average value of 60 naut. mi/deg 

1t T =+ 120 ' 

X =0* —- 0t - constant 15* = tor Case II and III 

= +30* ------ sketched 



FIGURE 4. 2 APPROXIMATE FOOTPRINT FOR THE STUDY VEHICLE 
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ALTITU0E, < ~KJ* FEET 

FIGURE 4. 4 TYPICAL NORMAL FORCE CHARACTERISTICS 
FOR HYPERSONIC GLIDER 
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SECTION V 

SOLUTIONS OF THE EQUATIONS OF MOTION 

This section deals with solutions of the complete equations of motion as de¬ 
rived in Section III. Solutions described here represent: (1) three basic cases 
(Cases I, II, and III) depending upon the initial heading of the vehicle and the bank 
angle maintained during re-entry, and (2) solutions which represent perturbations 
in control input from the third basic case. Section 5. 1 defines the general conditions 
used for solution of the equations of motion including vehicle, geophysical, aero¬ 
dynamic, and control system characteristics. The solutions, themselves, includ¬ 
ing perturbations, are presented in Section 5.2. Case III involves a roll to a 45® 
bank angle early during the re-entry, and special attention is given to this maneuver 
in Section 5.3. 

5 .1 CONDITIONS FOR SOLUTION 

The example vehicle and the specific data and equations used in computing 
the re-entry solution presented in this report are described below. 

5.1.1 Vehicle Configuration 

A typical configuration for a maneuvering re-entry glider was assumed for 
the study vehicle. The configuration is depicted in Figure 5.1, and the following 
physical characteristics were assumed for use in the equations of motion. 

= 23.6 psf IlS 1 * = .94 = 0 

I'~ J l = -.66 =-.74 

J x TZ 

The principal axes, the body axes, 'and the reference line for measuring angle of 
attack, etc. are all assumed to be coincident. 

5.1.2 Characteristics and Range of Solutions 

Typical re-entry solutions were computed which included the following 
general characteristics: 

a. Flight path skewed to the equator, 

b. Skipping flight path, 

c. Turn to produce lateral range, 

d. Rapid maneuvers to produce large rotations about all three vehicle axes. 

In order to appraise the validity of simplification in the equations, it is necessary 
to study their effect in the most general type of motion and not ir. a special case. 

Case III was just such a general solution, and provided a good medium for investigat¬ 
ing simplifications. 

The following initial conditions, typical for a near-orbital hypersonic 
glider starting re-entry, were common to the three basic cases. 

4 = 300,000 ft 

r. = -i* 

y = 24, 000 ft/sec 

a. = 15® 

The initial conditions for the three cases differed in the following characteristics. 

Case I: =0®, P a =0®; Case II and III: 0^ = 120®, P =30® 

The method of computing all the required initial variables is given in Appendix D. 
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The solutions were programmed so that they simulated tasks that might be 
asked of the pilot. This was done by commanding the control system to hold the 
following angle of attack ), sideslip ), and bank angle (0 ) time histories. 

Case I and II: £ = 15*, ^ = 0*, 0 = 0* 

Case III: <2 = 15*, 0 = 0 # , 0 = 0* to t = 185 sec 

& - 15*. /6 = 0*. 0 = 45* for t > 185 sec 

InCase III at "t - 185 sec a roll to 0 - 45* is commanded. This time 

corresponds to the bottom of the first skip, or the first occurrence of level flight 
in the re-entry. The roll maneuver produces large motions about all three vehicle 
axes (due to inertia coupling), and thus provides a good medium for studying the 
short-period or rapid modes of motion. The commanded angle of attack (dfe = 15*) 
was selected to correspond roughly to that for (L/D) max . Thus the three basic 
cases cover a wide range of re-entry tasks, and at the same time compare the 
effects of maneuvering, turning, skipping, and varying the flight path direction. 

The solutions were terminated according to the criterion A < 100,000 ft. 
The Mach number at this point was roughly M = 2. 5, so that the vehicle was 
now flying at a conventional altitude and supersonic speed. This point was also 
chosen as the terminating point because it is about the beginning of the landing 
approach. 

5.1.3 Equations of Motion 

The solutions appearing in this section (V) were obtained using the complete 
equations of motion as derived in Section III and shown in fully expanded algebraic 
form in Appendix A. The specific equations solved by the IBM 704 computer are 
all the equations listed in Appendix A, which are numbered to correspond with 
tKe'same equation in matrix form in Section III. The computer program used 
the matrix form of the equations directly, and the following listed equations and 
associated scalar and matrix definitions form the complete set. 

Equations: 3.37 - 3.51, 3.55 - 3.57, 3.55s, 3.56s, 

3.59 - 3.64, 3.68 - 3.71, 3.73 - 3.81 

Though it was intended to include wind and gust data in the equations, this 
was not actually done for lack of time. Also, the applied moments given by 
Equation 4. 39 were inserted as [My] in Equation 3.78. A three-axis flight con¬ 
trol system was simulated, as defined by Equations 4. 39, to provide vehicle path 
and attitude control. This flight control system, which because of the numbers 
used might be thought of as a combination of pilot and automatic control system, 
is described specifically in Section 5. 1.6. 

Numerical integration was performed using a fourth order Runge-Kutta 
method. The time interval or step size was 5 sec, except for 21 one-second 
steps which were used at t" * 185 sec when the 45* right bank was commanded 
by the controller. This numerical integration method and other pertinent digital 
techniques are described in Appendix D. 
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For a ll but three time history solutions which were stopped a.tt = 500 sec, 
the digital program stopped computing on meeting the condition A 100,000 ft. 

5. 1.4 Geophysical Characteristics 

The geophysical characteristics of the earth used for solution of the 
equations of motion in this chapter have been described in Section 2. 2 and 2.3. 

For convenience, however, the specific characteristics are included here also. 

The earth is assumed to be spherical, to be rotating at a constant 
angular velocity, and to have an inverse-square gravitational field. The following 
constants are used: 

1. Radius of the earth, = 20,860,000 ft 

2. Angular velocity of rotation, -XI = 0.004178 deg/sec 

3. Acceleration of gravity at 

the earth's surface, <3 0 = 32. 17 ft/sec2 


The atmospheric properties used in the solutions are the density and 
temperature. The density is given by the equation: 

/> •/>„ « : JK7 ^ T 

in which X = 0 , 0 * ^ < 25, 000 ft 

= 0.2 £ - 5000 , 25, 000 ft < £ 65, 000 ft 

= 8000 , 65,000 ft «< < 


(5. 1) 


and the sea level density is 

= 0.002378 slugs/ft 3 

The atmospheric temperature variation with altitude is programmed as a table 
which is linearly interpolated by the computer: 


^ ~ Feet 

7*~ °/f 

0 

36,500 

82,000 
156,000 
175, 000 
250,000 

518.688 

389-988 

389-988 

508.788 

508.788 

354.35 


The wind velocity model was programmed, but it was not used. No 
gust velocity model was programmed. 

5.1.5 Aerodynamic Characteristics 

General aerodynamic characteristics of a re-entry vehicle are described 
in Section IV. It is obvious that a complete aerodynamic representation would be 
very complex; consequently, to facilitate solution of the equations of motion as 
presented in this section (V), a much simplified representation is used. 

Trimmed lift and drag coefficients were used to represent the longitudinal 
forces. This expedient simplifies the force equations, eliminates the need for 
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considering elevator deflections as such, and gives good over-all values for 
lift and drag. Though errors result in some of the force derivatives, notably 
and , these errors are not significant to the results of the study. 

The representation is given, then, in terms of # and f body axis trimmed force 
coefficients, C Xjm and Cj t . Expressions for these coefficients are as follows: 

*■ C Xy . = 0,025 

= -0.0325 + (0.401 + 2.063-jjj. ) & 

+ (1. 166 - 2.298 ^ )\<t j A 


(5.3) 

(5.4) 


in which Aj is the Mach number, and <k is the angle of attack in radians. 
These expressions are plotted Figure 5.2a for dr in degrees. Aerodynamic 
forces are obtained by inserting Equation 5. 3 and 5.4 in the matrix equation: 

N * [<■] (S.5) 


in which 



(9.6) 


The aerodynamic moment coefficient matrix [C A was set to sero, and all 
moments were handled through [>V r ] as described in the next section. 


To provide a better visualization of the trimmed longitudinal aerody¬ 
namic representation, the trimmed lift and drag coefficients, c L r and C fir , 
and the trimmed lift-drag ratio, (L/0) r , have been obtained and are plotted 
in Figures 5.2b and 5.2c respectively. The matrix equation 


.<]•**[£] 


(5.7) 


defines the transformation from body axes to wind axes when the sideslip, , 
is zero. 


5, 1.6 Control System - Concept and Characteristics 

In order to compute solutions of the equations of motion which are 
similar to what actual re-entry flights might be, some method for controlling 
vehicle attitude and path had be be devised. The normal process would have 
been to incorporate the complete aerodynamic characteristics for the vehicle 
into the equations, and then design and incorporate into the equations a flight 
control system for the vehicle which would control the aerodynamic surfaces 
(elevator, ailerons, and rudders) and reaction controls so as to provide the 
vehicle with a desirable response to pilot command inputs. The flight control 
system would have logically been an adaptive one, which would have provided 
a sensibly constant dynamic response throughout the re-entry. For this study, 
an expedient short cut was taken. A desired response in (k »/0 , and 0 was 
hypothesized, and then only the moments necessary to produce this response 
were applied to the vehicle. The selected control scheme, involving six 
feedback loops, was one that could be easily incorporated into the equations, 
had physically significant parameters, and simulated roughly what might be 
realized from a good adaptive control system. 
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The three-axis controller, described by Equations 4.39 - 4.42 in 
Section 4. 5, provides several of the conventional stability derivatives, namely 
» £**4 , c * f . and , and the derivative is similar to the 

conventional Cji ^ . The C\g derivative is not one that is inherently present, 
and would only be formed as an element of the flight control system. In 
addition to its primary function of attitude stabilization and control, the 
three-axis controller provides path control capability for the vehicle. Thus, 
the & controller allows variation of the lift-drag ratio with the consequent 
variations of range, maximum heating rate, total heat input, deceleration, 
etc. The 0 controller can be used to rotate the total lift vector about the 
X body axis, thus providing capability for lateral or cross range movement 
of the vehicle. Most solutions presented in this section (V) makes use of the 
df controller to maintain a nearly constant 15* angle of attack corresponding 
roughly to (L/D) max for the programmed C% and 6 ^ . Case III solutions 
include a banking right turn which is accomplished by rolling the vehicle to 
0 = +45* by commanding a step in 0 C , and perturbations from the basic 
Case III solution involve step commands in d! fi , , and 0 C . 


The three-axis controller is described by the following matrix 



L 


(0-*e) * jf 


M 

» 



/V 


(*-/*•)+/* 

In this equation, the various 4 

's are gain constants; 


(5.8) 


a c , and ^ are 

command inputs to the controller; the three variables ( 0 , dt , and /g ) and 
their rates are available in the computing program. The gain constants were 
selected to produce desirable natural frequencies and damping ratios for the 
equivalent linearized equations of motion. In this study, the natural frequencies 
and damping ratios were made respectively equal for all three control modes. 


In order to shorten computation time for the solutions, the highest 
frequencies present in the equations of motion (those due to the controller) 
were set at O*^ =0.3 rad/sec, giving a period of about 21 seconds. This 
allowed an integration step size or time interval of 5 seconds. The frequency 
( = 0. 3 rad/sec) is approximately that of the longitudinal "short" period 

mode given by Etkin (in Reference 1.6, Figure 8 ) for an altitude of 250,000 ft. 
The controller gains were selected to give a damping ratio ? = .7. The 
gain constant follows. 

*0 9 TT = -0.09, ** * TT = -0.42 

= 4?* = -0.09. = -0.42 (5.9) 

= Xi =+0.09, =+0.42 

5. 2 RE-ENTRY SOLUTIONS 

The intent of this project was to develop equations of motion for use 
with a simulator and then to compute re-entries to check the equations, to 
determine the effect of small changes in flight conditions (perturbations), 
and to determine what mathematical or geophysical simplifications, if any, 
could be made in the equations of motion (simplifications). Because of 
unavoidable delays in checking out the computer program, the numberof 
re-entry solutions for perturbations and simplifications was extensively 
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curtailed. Nevertheless, this section (V) and Section VI clearly demonstrate 
the capability of obtaining excellent re-entry solutions with the complete 
equations of motion, and the two sections(V and VI) present a number of solu¬ 
tions indicating effects of the more important perturbations (from the piloting 
task point of view) and simplifications. 

Section 5.2.1 intorduces the three basic cases (I, II, and III) chosen for 
analysis and describes perturbations made on Case III. Section 5.2. 2 presents 
the solutions for Cases I, II, and III, and Section 5.2.3 presents the perturba¬ 
tions of Case III. 

5.2.1 Introduction 


Three basic solutions were computed during the course of this program: 
Case I, a re-entry in the equatorial plane, was chosen as the simplest to chpck 
and analyze and to provide a comparison with Case II, a re-entry whose initial 
heading is roughly along the Atlantic Missile Range. Cases I and II have identical 
initial conditions except for the heading angle. Case III may be compared with 
Case II to show the effect of a banking turn upon the re-entry trajectory. It is 
probably a more meaningful re-entry solution insofar as pilot and control system 
tasks are concerned because of the banking turn and because deceleration, heating 
rate and dynamic pressure are greater. In addition, Case III illustrates inertia 
coupling during the bank angle transient response and the lateral or cross range 
capability of the winged re-entry glider. 

As far as initial conditions are concerned, the three basic cases are nearly 
identical. For all solutions, the geophysical representation of the earth, the 
controller parameters and the vehicle's aerodynamic representation are identical. 
These are described in Section 5.1. Initial angular velocity of the vehicle in 
inertial space is zero. Angle of attack is initially set at +15*; the control system 
attempts to maintain this value throughout the re-entry except, of course, during 
Qc perturbation solutions. Large deviations from the commanded , and 

/S occurred only during the roll maneuver, caused by inertia coupling. The 
initial sideslip is zero; therefore, the condition & = holds and yields an 

initial pitch angle of +14*. Initial sideslip is zero; therefore, the Euler angle 
& is initially equal to the heading angle V r . In every solution the initial 
bank or roll angle is zero. Since it makes no difference to the solutions, the 
initial longitude in every case is zero. These are the initial conditions which 
are identical in every solution. 

The differences among the three basic cases are best shown in tabular 
form with some initial trajectory parameters -altitude, inertial velocity and 
flight path angle— included for completeness. 
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RE-ENTRY 

INITIAL CONDITIONS 

.— 

CONTROL INPUTS 

CLASSIFICATION 

A 

ft 

ft/sec 

deg 

deg 

* 

deg 

DURING RE-ENTRY 

Case I 

"Equatorial" 

300,000 

24,000 

- 1 

0 

+90 

&C = +15* 

4 = 

0c = 0 

Case II 

"Atlantic Missile 
Range" 

300,000 

24,000 

-1 

+ 30 

+ 120 

= +15* 

/Sc = 0* 

0c =0* 

Case III 

"Atlantic Missile 
Range" with Banked 
Turn 

300,000 

1 

24,000 

-1 

+30 | 

! 

+ 120 

^6 = +15* 

0c = 0‘ 

For t<■ 185 sec 

0c = 0* 

1 For 185 sec 

0t = +45* 

(Step input in 0 C ) 


The table shows that the basic differences in the three cases are initial values of 
latitude and heading angle, and a banking maneuver during the solution. 


Since the computer program input routine requires the initial values of those 
twelve variables whose rates are integrated during solution ( U , sir , ur ; jO , A , /* i 
, 0 ; A. , , X ), specifying parameters such as V , 2f , ft -» OL , etc. is not 

sufficient. Appendix D shows a sample calculation for the required values of the 
computer variables for Case III. In the general case, the computation of initial 
conditions is not so easy as shown for Case Ill, a case which had no initial angular 


rates, bank angle or sideslip. 


A re-entry trajectory is fraught with limits pertaining to maximum heat 
input, maximum heating rate, maximum dynamic pressure, maximum decelera¬ 
tions, maximum lift, etc. These limits define a re-entry corridor which is 
usually narrow; therefore, it is of interest to present trajectory solutions in a form 
which can be related to the re-entry corridor. To this end, the data is presented 
in graphical and tabular forms; each is necessary to illustrate certain aspects of 
the re-entry solutions. Especially important from the re-entry corridor point of 
view are the altitude (A. ) versus relative velocity ( V £ ) plots,. In general, the 
re-entry solutions presented here stay within reasonable limits for a gliding 
re-entry vehicle. This is not surprising since the vehicle has a high lift-drag 
ratio ( ** 2.4), the initial inertial velocity (24, 000 ft/sec) is lower than the 
circular orbital velocity (25,721 ft/sec at the same altitude with the same sea 
level gravity), and the initial flight path angle is small (-1*). 


Small changes in the commanded attitude provide perturbations of the 
basic trajectory. Case III. These solutions are designed to show how carefully 
a given attitude must be held to result in desired range, cross or lateral range, 
heat input, deceleration, etc. The results are presented in tabular form; 
however, where the perturbations produce significant differences in certain 
variables, a graphical presentation is also used. 
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5.2.2 Presentation and Comparison of Three Basic Cases 

The major trajectory parameters for the three basic solutions are presented 
as follows: 


Case I Figure 5. 3 

Case II Figure 5. 4 

Case III Figure 5.5 

Each figure presents time histories of altitude £ , relative velocity Vg , and 
flight path angle # , and a plot of vehicle altitude versus its relative velocity. 

Since the coordinates of the n and Vg plot are closely associated with potential 
and kinetic energy respectively, this plot is useful from the energy management 
point of view. In fact, altitude versus relative velocity would be a desirable display 
for the pilot of a re-entry vehicle especially if the instrument had a re-entry cor¬ 
ridor overlay. The important Case III presented in Figure 5.5 also includes a 
map-like plot of latitude ( -j) ) versus longitude (A ) which is the vehicle's track 
over the surface of the earth during re-entry. 

Comparison of the three basic cases is made using Figure 5.6 and 
Table 5. 1. The figure contains the same data as in Figures 5. 3 through 5.5, but 
each plot contains three curves, one for each case. The table presents trajectory 
data for three distinct conditions during the re-entry: 

1. The bottom of the first skip; that is, the first point on the 
altitude versus time plot at which altitude reaches a mini¬ 
mum point. 

2. The "top" of the first skip; that is, the first point at 
which altitude reaches a maximum. 

3. The point at which the vehicle reaches an altitude of 
100,000 ft. This is the condition for stopping the 
particular solutions. 

None of these conditions is met exactly because of the discrete-time nature of the 
solutions, a five second time interval being used for numerical integration. 

Certain key variables are tabulated which show the variation from the exact desired 
conditions. For conditions (1) and (2) above, the flight path angle ? is the tabu¬ 
lated key variable. Whenever altitude is at a maximum or minimum, the flight 
path angle should be exactly zero. Similarly, altitude is the tabulated key variable 
for condition (3) for which A should be always exactly 100,000 ft. In view of these 
considerations, perfect comparisons using the tabulated data cannot be made; 
however, as the indicators ( A and # ) show, the conditions are met closely 
enough so that valid comparisons can be made. 

The values of range and stagnation-point heat-transfer rate tabulated in 
Table 5. 1 were hand calculated using computer print-out results. Range, in 
nautical miles, was calculated from the initial and final values of latitude and 
longitude using the equation 

&AM6E - Ccs m/ £s/>! *" 00 % 7% cos % (5.10) 

where is in nautical miles. Care must be exercised in evaluating the arc¬ 
cosine. The stagnation point heat transfer rate was calculated using the equation 
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Q=Z2x/o- 9 \fir Vj ~sec* (5.11) 

in which /o and I/*, are expressed in slugs/ft 3 and ft/sec respectively. /? is the nose 
radius of the body in ft, and /? = . 5 ft was assumed for the nose of the study vehicle. 

5.2.3 Perturbations - Case III 


Table 5.2 presents the results of introducing small attitude perturbations in 
Case III at 185 seconds* of the time history. These perturbations in the attitude 
command were held throughout the rest of the solution. To focus attention upon the 
linearity (or non-linearity, as the case may be) of the perturbational effects on 
trajectory parameters, the changes in the parameters rather than the parameters 
themselves are tabulated in Table 5.2. Since the commanded changes in <2 , yg , 
and 0 were small (±1* from conditions of Case III), the differences in the trajec¬ 
tory are generally too small to be presented in time history plots. This is another 
reason for presenting tabulated differences between the perturbation solution and the 
basic Case III. 

The perturbation solutions are six in number: three, those with increased 
O', £ , and 0 , were stopped when the altitude became 100,000 ft; and three, those 
with decreased Ot ,^ , and <0 , were stopped after 500 seconds of the time history 
had elapsed. Inadvertently, the latter three solutions were obtained as perturba¬ 
tions of Simplification Solution III-S1 (III-S1 is described in detail in Section VI). 
Since this simplification solution closely approximates the complete solution, it is 
possible to consider all six perturbation solutions as perturbations of the basic 
Case III. Nevertheless, in Table 5.2 the differences listed for Solutions III-P4, 
III-P5, and III-P6 were obtained by using Solution III-S1. This is a third reason 
for presenting tabular differences. 

One result has been presented graphically in Figure 5. 7; it is an altitude 
time history comparison of Solutions III, III-P1, and III-P4. This figure shows 
the angle of attack perturbations for the first 700 seconds of time history. The 
0L - 14* case (Solution III-P4) ends at 500 seconds, the program stop for that 
solution. 

So far as conditions at solution termination are concerned, the 0 and 
perturbations have a larger effect than the A perturbation. This is probably 
because near Ot for (L/D)max. small changes in A have no effect on L/D. 

Table 5.2 contains data which compares the perturbation solutions III-P1, 

III-P2, and III-P3 with III when A 100,000 ft, at which condition the small 
attitude perturbations have had sufficient time to produce noticeable changes in 
the trajectory. 

A principal reason for computing the perturbation solutions was to obtain 
a measure of the trajectory controllability, and hence to provide a yardstick for 
measuring the significance of errors introduced by simplfications in the equations 
of motion. 


* The same time that the banking turn is introduced during each solution. 



NAVTRADEVCEN 801A 


5.3 ROLL MANEUVER AND TURN DURING RE-ENTRY 
5.3.1 Introduction 


The turn in the Case III solutions is important in this study because it 
allows evaluation of the maneuverability of the vehicle. The roll entry into the 
turn is especially important because, due to inertia coupling, dynamic responses 
in all three of the attitude modes of motion are excited, and this represents a very 
general type of short-period response. 

5.3.2 Case III - The Turning Re-Entry 

Trajectory parameters of the basic Case III solution are presented in 
Figure 5. 5 and Table 5. 1 in sufficient detail that discussion here can be brief. 

Figure 5.6 provides a graphical comparison of Cases II and III which differ only 
in the commanded bank angle. The one quantity of interest that is not numerically 
specified in Figure 5..6 is the lateral or cross range achieved by the turn of 
Case III with respect to the straight-in solution of Case II. Roughly, the cross range 
achieved is 1600 nautical miles, certainly a significant capability. This informa¬ 
tion was obtained from the map-like Figure 5.6c which compares the latitude versus 
longitude plots for the three basic cases. 

It is readily apparent that cross range obtained by rotation of the lift vector 
is costly from other points of view. The effective (vertical) lift-drag ratio is 
considerably smaller by some 30% for the 45* bank angle, resulting in increased 
deceleration, heating rate, and dynamic pressure; for as Figure 5.6b shows, at any 
given altitude (consequently, air density) after the roll angle input, the turning 
re-entry of Case III has a greater relative velocity than Case II. 

The roll input command ( = +45*) at t - 185 seconds produces a dynamic 

response of the vehicle in roll which is governed essentially by the 0 and 0 
derivatives of the three-axis controller described in Section 5.1.6. Equally 
significant is the transient in dr and /S caused by inertia coupling. Figure 5.8 
presents time histories of dr , /6 , and 0 for the period of time which includes the 
response transient. 

5.3.3 Perturbations 


From the piloting-task point of view, the perturbations of Case III, 
presented in Table 5. 2, shed some light on the care with which the vehicle's 
altitude must be controlled to produce a desired trajectory. Suppose that a given 
angle of attack must be held to achieve a desired cross range. Does a one degree 
change in a specific attitude angle seriously affect energy management? Will 
the pilot's instruments have to be accurate enough to indicate an angle within a 
few tenths of a degree or will less precise instruments suffice? How much 
steady-state error can be tolerated in an attitude control system before the 
system ceases to be adequate? A rough answer to such questions can be made 
based on the perturbation solutions. 
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CONDITION 

VARIABLE 

I 

Equatorial 

Re-Entry 

11 

Atlantic 
Miaaile Range 
Re-Entry 

III 

Turning 

Re-Entry 


A 

ft 

206,677 

207,337 

207,337 

Firat Minimum 

Kr 

ft/aec 

22,082 

22,385 

22,385 

Altitude Point, 

it 

deg 

-0.0034 

40.0597 

40.0597, 

"Bottom" of 

t 

aec 

185 

185 

185 

Firat Skip 

a 

BTU/ft 2 

228 

240 

240 


h* 

lb/£t^ 

91.2 

91.1 

91. 1 

Firat Maximum 

A 

ft 

276,236 

276,312 

252,278 

Altitude Point, 


ft/aec 

21,442 

21,442 

21,635 

"Top" of Firat 

It 

deg 

+0.0046 

-0.0213 

40.0083 

Skip 

t 

aec 

390 

390 

366 


f;* 

lb/ft 2 

4.5 

4.5 

12.7 


A 

ft 

100, 184 

100,148 

100,032 


* 

ft/aec 

2,033 

2,005 

2,464 


it 

deg 

-3.86 

-3.65 

-4.03 

End of Solution 

A 

deg 

4122.53 

4113.66 

449. 17 

A 2* 100,000 ft 


deg 

0 

-33.66 

-31.11 


Pr 

deg 

490.00 

4117.00 

4 350.29 


t 

aec 

3.035 

3,065 

2,146 


Range 

naut. 

mile 

7,330 

7,450 

4,607 


* f i the incompreaaibla dynamic preaaure. 


TABLE 5. 1 COMPARISON OF THREE BASIC CASES 



TABLE 5. 2 EFFECT OF ATTITUDE PERTURBATIONS ON 
TRAJECTORY PARAMETERS 
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FIGURE 5. 1 THREE VIEW OF TYPICAL RE-ENTRY VEHICLE 
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a. Body Axes Force Coefficients 



b. Lift and Drag Coefficients 


rf 

a 

r 

? 

I 



c. Lift-Drag Ratio 

FIGURE 5. 2 LONGITUDINAL AERODYNAMICS - 
TRIMMED CONDITION 
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FIGURE 5. 4 ATLANTIC MISSILE RANGE RE-ENTRY SOLUTION 
CASE II 















NAVTRADEVCEN 801A 



FIGURE 5. 5 TURNING RE-ENTRY SOLUTION - CASE III 















TIME, t ~ SECONDS 


FIGURE 5. 7 


ANGLE OF ATTACK PERTURBATIONS 
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FIGURE 5. 8 SHORT PERIOD FIGURE 5. 9 SHORT PERIOD 

TRANSIENT, ROLL ANGLE TRANSIENT, ANGLE OF ATTACK 
INPUT - CASE HI PERTURBATION 



FIGURE 5. 10 SHORT PERIOD FIGURE 5. II SHORT PERIOD 
TRANSIENT, BANK ANGLE TRANSIENT, SIDESLIP 

PERTURBATION PERTURBATION 
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SECTION VI 

SIMPLIFICATION OF EQUATIONS OF MOTION 

6.1 INTRODUCTION 

This section (VI) is important to the successful use of the equations of 
motion for a hypersonic re-entry vehicle in a simulator with a real-time digital 
computer. The very concept of a piloted simulator implies precision of a low 
order relative to the precision required for pre-programmed satellite launchings, 
etc. The simulated short period dynamic motions must simply satisfy the perti¬ 
nent senses of the pilot "flying" the simulator. Trajectory parameters which are 
observed over long periods of time by the simulator pilot must be accurate, but 
not so accurate as calculation of orbital ephemeris, etc. Consequently, the 
equations of motion were investigated for possible simplification to facilitate real¬ 
time computation and integration of the equations, while providing adequate simu¬ 
lation in all long and short period aspects of re-entry. 

The following sections (6.2 and 6.3) will respectively deal with assumptions 
made to simplify initial derivation of the equations, and subsequent mathematical 
and geophysical simplification of the derived equations. The latter simplifications 
will be treated much more thoroughly in the following discussion because they 
involve solutions of the equations of motion and because the simplifications used to 
derive the equations are discussed in Section III and Appendices B and C. 

6.2 SIMPLIFYING ASSUMPTIONS FOR DERIVATION OF EQUATIONS 

OF MOTION 

There are three basic simplifying assumptions which were made to facili¬ 
tate the initial derivation: one is geophysical in nature and the other two concern 
the vehicle configuration. It should be noted, of course, that there are other 
assumptions made during the derivation of Section III, but they result in neglecting 
relatively minor effects which are not at all well known or which do not pertain to the 
objectives of the derivation. The basic assumptions are listed: 

1. Earth oblateness is neglected, 

2. The vehicle's mass is constant, and 

3. Large thrust forces are not considered. 

6.2.1 Earth Oblateness 


A spherical earth is the basis for the derivation of equations in Section III. 
This is a reasonable simplification of the more nearly correct oblate spheroidal 
character of the earth which greatly facilitates the derivation and which results 
in a computation scheme economical of time. Parameters, needed to locate a 
vehicle relative to the spherical earth ( A , A , and/> ), comprise a set of pseudo- 
spherical* coordinates whose rates are integrated during the computation process. 
Aside from the purely mathematical convenience afforded, the spherical earth 
assumption appears to be adequate for use with a simulator when one considers 
some reasons for using the simulator. For example, familiarization with vehicle 
controls, instruments and handling qualities is not affected, crisis procedural 


* The description pseudo-spherical is used because the actual spherical 
coordinates would be the quantities + A , X and 90* - ^ respectively. 
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training is not affected, and vehicle dynamic response (as indicated by instrument 
panel displays) is not affected. These are among the most important considerations 
for good simulation. 

Flight time of a gliding re-entry vehicle can be substantial, especially if the 
vehicle is required to orbit the earth once, twice, or more often. Under these condi¬ 
tions it may be important that the pilot perform navigational tasks, thus requiring 
that consideration be given to navigation during pilot training. To obtain accurate 
computation of the trajectory or path of the vehicle for long duration flights, earth 
oblateness and related effects should be included in the equations; however, this 
addition is necessary only to obtain trajectories which are to match real trajec¬ 
tories. Navigational techniques used in a simulator would be identical with either 
earth model, although the paths realized would be different; therefore, there is no 
justification for the added complication of earth oblateness as far as navigation 
training is concerned. Another reason for neglecting oblateness is that the equa¬ 
tions are developed basically for re-entry which is a short duration, short range 
phase of flight relative to the orbital phase. 

The most difficult objection to overcome, with regard to neglecting earth 
oblateness, is the effect of oblateness on large-scale maneuvers. For example, 
if a control input had been computed to provide a specific resulting location for 
the real vehicle based on the oblate shape of the earth, to obtain the desired 
result with a simulator would require equivalent modelling of earth oblateness, 
all else being equal. Thus, where accurate simulation of the trajectory over a 
large range is required, earth oblateness will need to be included in the equation. 

The required equations are derived in Appendix C to fill this need. Nielsen, 

Goodwin, and Mersman (Reference C.4) give solutions from which a measure of 
the effect of oblateness can be obtained. 

6.2.2 Thrust and Variable Mass 


While earth oblateness is a geophysical aspect of the simplification 
problem, consideration of variable mass and thrust are strictly aspects of the 
vehicle itself. The equations of motion developed in Section III are capable of 
handling the near-orbital phase of a boost-orbit-re-entry flight as well as the 
re-entry phase, but main emphasis for the derivation was placed on the re-entry 
phase. Since for re-entry only the glide vehicle is involved, the mass is very 
nearly constant (the small loss being due only to reaction controls). It is likely 
that the glider may use a landing engine, but this has not been specifically con¬ 
sidered in the derivation since methods for handling such engines are well 
known. Consequently, the vehicle has a constant mass and no thrust, assump¬ 
tions founded in the ground rules of the study. 

If the desired configuration is composed of booster and glider, variable 
mass and large thrust forces certainly must be considered. The allied prob¬ 
lems of variable moments of inertia, inclusion of products of inertia which vary, 
and moments due to the thrust also require attention. These problems are 
treated in Appendix B. 

6.2.3 Less Significant Assumptions 

It should be mentioned for completeness that other, less significant, 
assumptions are made for the derivation of Section III. Most of these are made 
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because of lack of knowledge of a particular phenomenon, a situation apparent in 
mathematical representation of the geophysical environment, especially winds and 
gusts. Other assumptions, such as considering a non-rotating reference frame 
fixed at the earth's center as an inertialframe for the motion of the vehicle rela¬ 
tive to the earth, are made because they have-completely negligible effect on 
solutions of the equations of motion. 

6.3 EFFECT OF SIMPLIFICATIONS IN THE EQUATIONS OF MOTION 

Two types of simplifications are dealt with in this section. Although both 
types are accomplished in essentially the same manner, by adjusting the computer 
program, they are basically different in effect. The simplifications are: 

1. geophysical, or 

2. mathematical 

in nature. Presumably, simplifications of complete six-degree-of-freedom equa¬ 
tions of motion can be made in many ways. The simplifications described in this 
section were made: first, by eliminating terms in certain equations; secondly, by 
forcing certain constants to zero; and thirdly, by making certain variable quanti¬ 
ties constant. Linearization, a widely used technique for simplification of mathe¬ 
matical equations, was not explicitly attempted in this study although it could for 
special cases be advantageous. An example of the sort of linearization which 
might be used is the replacing of sines and cosines by their small angle approxi¬ 
mations . 

6.3.1 Mathematical Simplifications 

There are three basic mathematical simplifications described in this 
section. They involve two successive approximations of a lengthly equation, an 
approximation for the vehicle's angular velocity with respect to the earth, and 
simplification of the equation for the rate of change of the vehicle's orientation 
angles. It is not implied that these exhaust the possible mathematical simpli¬ 
fications. The simplifications examined were those thought to be most desirable 
for facilitating real-time digital computation. 

The most lengthly equation derived in Section III is that for the rate of 
change of the vehicle's velocity with respect to the air (see Equation 3. 73 and 
Equation A. 73). Equation 3.20£ with wind and gust terms omitted, is used to 
explain the first two simplifications: 

(V*) s * +Z* *V tr ( 6 . 1 ) 

Solution III-SI shows the effect of neglecting the (Vtr)g term in 
Equation 6. 1. Trajectory data is presented in Table 6.1 so that parameters may 
be compared with Case HI. Table 6.2 must be used to determine the effects of 
the simplification upon the short period dynamic response because the differences 
between III-S1 and HI are too small to illustrate graphically. 

A further and final simplification which can be provided for Equation 6. 1 
is to eliminate the term in addition to the term. This simplification 

demonstrated by Solution HI-S2, is presented in Tables 6.1 and 6.2 and Figure 6.1 
Comparison of HI, IH-S1, and IU-S2 shows that all three solutions are nearly 
identical as far as trajectory parameters are concerned. Figure 6.1 shows that in 

*see page 47 
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t 


attitude dynamics III-S2 is considerably different from III. Solution III—S1 was not 
plotted on Figure 6. 1 because no difference between it and III could be shown. 

The results of these first two mathematical simplifications inchoate that 
they are both good for trajectory purposes; but, whereas neglecting (7 er t>^ provides 
a good approximation for.the short period response, neglecting uf g x is 
definitely not allowable. On this basis, we can conclude that lVVerVl an< * 
may be neglected in Equation 3.73. L J 

Equation 3. 15 gives the vehicle's angular velocity relative to the earth. 

It is repeated here: 

UJ e * Zo - si (6. 2) 

Solution III-S5 indicates that neglecting -H- in Equation 6.2 provides totally 
unacceptable results in both the trajectory and short period responses (Figure 6.2 
and Tables 6. 1 and 6.2). Neglecting the earth's angular velocity in Equation 6.2 
is an example of an inconsistent geophysical simplification; it essentially removes 
earth rotation from the vehicle's rotational or angular equations while not 
removing it from translational or position equations. There are two specific 
effects of this_situation which bear mentioning. First, the odg-Eo simplification, 
like that of _/2 =0, produces a negligible effect on short period response 
(Table 6.2). Secondly, for the solution termination condition ( A =** 100,000 ft), 
cii £ » co produces effects of opposite sense from those of 71 - 0 in altitude A , 
relative velocity V £ , flight path angle # , as well as the vehicle's track 

(compare Figures 6. 2a and 6. 2b with 6.4a and 6. 4b respectively for details). 

Since the Zo £ =-co simplification affects the vehicle's orientation directly through 
Equation 3.55 and indirectly through Equation 3.73 and the control system, it 
is not difficult to see why for III-S5 is so much different from f for III 
(Figure 6.2a). 

The final mathematical simplification involves the Euler angle calculation 
given vectorially by Equation 3. 52: __ 

^ o- 6 + jj * co - Si - (%+■ -j>) (6.3) 

Solution III-S6 was obtained for the condition that 

f* + 6 -h J - (6.4) 

Table 6.2 shows that there is very little difference between the short period 
responses of III-S6 and III. This situation, like that of the Lo t =£i simplification, 
is a result of the fact that the short period response is primarily dependent upon 
the vehicle's physical characteristics and the three-axis controller. Table 6. 1 
indicates that the trajectory is very poorly computed with the condition of 
Equation 6.4. It is a completely worthiess_simplification, so poor in fact that 
the solution has not been plotted. If only SL had been neglected in Equation 6. 3, 
a more valid simplification might have been obtained. 

6.3.2 Geophysical Simplifications 

Two geophysical simplifications were investigated-constant gravitational 
acceleration and non-rotating earth. It would seem that these are rather drastic 
simplifications of the vehicle's complex geophysical environment; however, they 
are the basic assumptions made for analytical solution of re-entry trajectory 
equations. How reasonable are these geophysical simplifications from the 
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simulator point of view? The tabular and graphical results presented in this 
section should satisfactorily answer the question. 

An inverse-square gravitational field is programmed with the equations of 
motion (see Equation 3. 81). At an altitude of 300, 000 ft, the initial altitude for 
all re-entry solutions, the acceleration of gravity is 

& =32.17 H P: 86 x = 31.2643 ft/sec 2 (6.5) 

(21.16 x 106)2 

It is this value of gravitational acceleration that is used in Solution III-S3 so 
that III-S3 will have initial conditions matching all other solutions. Solution III-S3 
is actually a simplification of Solution III-S2 since ( V a J_ = (V) B in both; thus the 
two are compared for trajectory parameters in Figure o. 3 and may also be com¬ 
pared in Table 6. 1. Constant gravity has too small an effect on the vehicle's 
short period response to be plotted; in fact, the little difference that might exist 
between Solutions III-S2 and III-S3 is not perceptible in Table 6.2 which contains 
data to just two decimal places. 

The graphical data shows that a constant g equal to the value at 300, 000 ft 
results in a small drift of the trajectory parameters from those of III-S2. The 
maximum difference between the two solutions could be lessened if some mean 
value of gravitational acceleration were used. Such a mean value can probably be 
determined for any desired re-entry simulation. However, in a digital computer 
program the calculation saved hardly seems worthwhile. 

The motions relative to the earth are described in terms of the vehicle's 
inertial linear and angular velocities and the rotation of the earth. In many 
situations it is not necessary to consider earth rotation; however, this conveni¬ 
ence cannot be justified in the equations of motion considered here because the 
duration of each solution is long enough for the earth to rotate through several 
degrees (corresponding to several hundred miles on the earth's surface) and 
because the vehicle's velocity relative to the air would be different, resulting 
in different aerodynamic forces and moments. 

To be able to assess the effects of a non-rotating earth, Solution III-S4 
was obtained. Trajectory parameters are plotted in the two parts of Figure 6.4 
along with the corresponding data for basic Case III solution. It is important to 
note (Figure 6.4b) that earth rotation cannot be accounted for just by a correction 
in longitude ( X ). A numerical comparison of the trajectory can be made by 
using Table 6.1. Like the constant g solution, Solution III-S4 shows that effects 
of the =0 simplification on short period dynamic response are completely 
negligible. Table 6.2 tabulates various short period dynamic response param¬ 
eters for Solutions III and I1I-S4. 

Of the two simplifications described in this section, constant g is by far 
the better. Both simplifications are equally good in short period dynamic 
response characteristics, but comparison of Figures 6.3 and 6.4 shows that 
constant g is better for trajectory computation. This is notan unexpected 
result since the variation of g from zero to 300,000 ft of altitude is small, 
less than 3%. 
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Geophysical simplifications are discussed in Section VII and recommenda¬ 
tions are made regarding whether the simplifications should be adopted or whether, 
in fact, more accurate geophysical representations are needed such as described 
in Appendix C. 






























(a) 

FIGURE 6.4 EFFECT OF EARTH ROTATION ON TRAJECTORY 
PARAMETERS 
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SECTION VII 
CONCLUSION 


7.1 INTRODUCTION 

The purpose of this concluding section is to summarize the requirements 
for simulating the re-entry and landing of hypervelocity aircraft as determined 
from the study described in this report. Although the study emphasized those 
aspects of the problem associated with re-entry, consideration was given to 
the requirements for simulating final approach and landing. In addition, informa¬ 
tion was obtained relative to boost and en-route cruise or orbital flight. Basically, 
this study concerned only the requirements for the equations of motion, but the 
rest of the simulation had to be considered in order to perform this task. 

Special importance was attached to the required 1 piloting task, since the accuracy 
required in the equations of motion is directly a function of the combined response 
of the pilot and the aircraft-of how tight the loop is closed. 

The development in this section for the recommended equations follows 
the outline of the report: geophysics, axis systems and coordinates, aerodynamics 
and control, and numerical methods and allowable simplifications. The accuracy 
required depends on the piloting task and the purpose of the simulator. Is it 
enough to simulate just the trajectory (energy management) problem or just 
attitude dynamics, or is complete simulation of the control problem required? 

Is accurate navigation important? Or is it enough to let the pilot just re-enter 
and land the aircraft anywhere? In this section an attempt has been made to indi¬ 
cate what is required in the equations for each of these types of simulation. 

Computational accuracy is an important factor in solving the equations. 

No problem was encountered with obtaining sufficient accuracy in this study using 
the recommended equations and an IBM 704 to solve them. However, if the 
accuracy of the computer used in a specific simulator (such as UDOFTT) is less 
than that of the usual general-purpose digital computer (i.e., IBM 704, etc.), 
then ways of improving the accuracy of the computation are described, especially 
through the use of special coordinates for integrating the equations. 

Although considerable effort has been spent in this study trying to develop 
the best equations for the task, there are many areas where further study is 
obviously needed. Some of these are: numerical integration techniques, coor¬ 
dinates for improving the inherent accuracy, simplifications in the kinematics 
and dynamics, representation and effects of aerodynamic forces and moments, 
stability and control requirements, and expected dynamic characteristics par¬ 
ticularly at high angles of attack. 

7.2 GEOPHYSICAL DATA 

Recommendations for geophysical representation of the earth are made 
in four parts: 

1. Earth motions, 

2. Earth shape, 

3. Gravitational field, and 

4. Atmosphere. 
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If a detailed specification of a simulation task is presented, specific recommenda¬ 
tions can be made. Here the discussion is based upon the premise of complete 
simulation with additional comments being made when necessary to describe less 
accurate simulation. 

7.2.1 Earth Motions 

As mentioned in Section 2. 2. 1, the earth has numerous motions in 
inertial space-the most important being translation in its orbit and rotation about 
its axis. Mathematical justification has been given in Section 2.2.1 and 3.2.1 
for neglecting the earth's translation; therefore, the rotation remains as an impor¬ 
tant motion. Results presented in Section VI, especially in Figure 6.4, indicate 
that accurate trajectory simulation requires the presence of earth rotation. 

It is recommended that constant earth rotation be included for proper 
simulation of re-entry flight under all conditions. 

7.2.2 Earth Shape 


There are many possible approximations for the shape of the earth; notably, 
flat, spherical and oblate. For re-entry simulation the flat earth is ruled out 
entirely (with one qualification; see Section 7.5.2). A spherical earth, used for 
obtaining the solutions of Sections V and VI, is adequate for short range re-entries 
(e.g. , small L./ D); but, in fact, the re-entry solutions obtained during this 
program have ranges too long for the spherical approximation. An oblate 
representation of earth shape should be considered the minimum for re-entry 
simulation though the full representation needed for predicting a satellite ephemeris 
is certainly not needed. 

A first-order approximation of the oblateness effect, as described in 
Appendix C, is recommended. This extension, of equations developed in Section III, 
is easy to incorporate and will result in significant improvement in a general simu¬ 
lation of re-entry. It provides not only oblate shape for the earth, but also the 
corresponding oblate shape of the atmosphere. 

7.2.3 Gravitational Field 

It has been shown in Reference 1.6 and elsewhere that an inverse-square 
gravitational force field is necessary for the accurate description of the skipping 
(phugoid) motions of a re-entry vehicle. Except for the constant gravitational 
acceleration solution, the solutions of Sections V and VI were obtained using 
the inverse-square approximation. It appears that nothing less accurate than 
this should be used with the equations of motion. Two comments can be made 
concerning this statement, however: 

1. The constant gravity solution (see Figure 6.3) is surprisingly 
good and it could have been better had some mean value of the 
constant acceleration been used instead of the value for 

A = 300, 000 ft. 

2. For simulation of long-range re-entries or flights involving one 

or two orbits of the earth, the oblate gravitational field is required. 
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Since the discussion concerns a general simulation, it is recommended 
that the gravitational field be represented by the gradient of the potential function 
(Appendix C, Equation C. 12; zero and second harmonic terms only). 

7.2.4 Atmosphere 

Physical properties of the atmosphere are not well defined at this time. 

In fact, mathematical description of these properties under all conditions is a 
seemingly hopeless task. Thus, someone intending to make use of these proper¬ 
ties must turn to mean data such as contained in Reference 2.6, the ICAO 
Standard Atmosphere; although one recognizes that complete correspondence 
between this mean and the actual conditions in the atmosphere is unlikely. 

This study dealt with re-entries using atmospheric temperature, density 
and pressure models for altitudes below 300, 000 ft, a region in which the proper¬ 
ties are relatively well defined in the mean. However, the properties exhibit 
large percentage variations from the mean which depend upon season, time of 
day, large scale surface weather conditions, etc. Above 300,000 ft even the 
mean is not well defined and again variations are large depending on the conditions 
mentioned above, with added dependency on solar activity and other extra¬ 
terrestrial phenomena. 

Figures 2. 1 and 2.2 compare the models of temperature and density 
developed in this study with the ICAO Standard Atmosphere. The temperature 
model is an exact match; the density model is only an approximate match but, 
in view of the discussion in the preceding paragraph, it is certainly adequate. 

It is impossible to fit a single unchanging exponential function to the 
standard atmosphere density data if adequate density description is required for 
the entire altitude range of interest. Without going to undue complexity, very 
little improvement can be obtained over the density model given by Equation 


The temperature and density models described in Section 2.3.2 are 
therefore recommended. 

Atmospheric motions are, like density, temperature, etc., poorly 
defined except at low altitudes where mean data on winds and gusts is available. 
Weather, time of day, season, etc. affect intensity, direction and gradients of 
the motions, thus making mathematical description of winds and gusts difficult. 
This is unfortunate, for wind gradients have been found to significantly affect 
short period motions of aircraft. 

Typical averaged mean wind data, modelled for use with the equations 
of motion (Equation 2.14), can be useful, but good simulation should possibly 
have a wind velocity vector model which is more accurately defined relative 
to geographical location and altitude and which more accurately depicts jet 
stream location. 

No gust model was devised during this study, but certainly a model 
can be used to advantage to provide good simulation of short period response 
and handling problems especially at low altitudes during approach and landing. 
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Mathematical application of gust data to equations of motion can be handled in 
many ways; it is felt, however, that the method described in Section III (using 
[/*■•] ) is entirely adequate for general simulation of re-entry flight. 

Atmospheric Oblateness: All properties of the atmosphere are presented 
as functions of altitude so that if altitude is measured above the oblate earth and 
is used for determination of these properties, the atmosphere will be similarly 
oblate. Atmospheric oblateness would be, therefore, one of the geometric effects 
of earth oblateness. 

7.3 COORDINATES AND AXES FOR INTEGRATION 

The coordinates and axis systems derived in this report are believed to 
be best for the simulation problem. Some comments will be made, however, to 
justify this statement and to indicate that one system can actually be better than 
another as far as accuracy of numerical integration is concerned. 

The three primary sets of axes developed in the study are listed here: 

1. Earth axes. This is a frame whose origin is fixed at the earth center 
and which rota tes at angular velocity jQ. with the earth. The location 
of the vehicle's center of gravity with respect to the earth is given 
in terms of spherical coordinates of this frame ( Y* eg , X , -p ). 

Z. Earth-surface axes. This is a frame whose -y f3 plane is parallel 

to the local surface of the earth, whose is directed toward the 

center of the earth, and whose origin is positioned at and moves with 
the vehicle's center of gravity. Since the position of the earth- 
surface frame is given in terms of Y* e g , /l , and , it is an ideal 
frame with which to reference geophysical parameters of the earth, 
especially gravitational acceleration, winds, etc. 

3. Body axes. The origin of this frame is located at the vehicle's 

center of gravity. Its X axis is directed forward; its y axis, to 
the right; and its ^ axis down (as viewed by the pilot). This frame 
is of prime importance in the equations of motion because the vehicle's 
force and moment equations are integrated in this frame. The vehicle's 
orientation is given by the orientation of the body frame relative to the 
earth surface frame. These axes are the natural coordinates of the 
vehicle, for instrumentation, for accelerations, for defining physical 
characteristics of the vehicle; and they are also convenient for 
representing aerodynamic forces and moments particularly at hyper¬ 
sonic speeds. 

There are other reference frames used which are of less importance to the 
mathematical development of the equations of motion. Two of them, the wind and 
stability axis frames, are used in obtaining and handling aerodynamic data. These 
frames have been adequately described in the Symbols and Definitions at the front 
of this report. The so-called trajectory ( T ) axis frame was not specifically 
used in the study but is, nevertheless, defined by the vehicle's center of gravity 
and the tangential velocity vector ( V T ) which is directed along the X T axis. 

Since it is conveniently present in the development of the equations, it was consid¬ 
ered as a possible frame in which to integrate the force and moment equations 
rather than the body axis frame. 
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It is important that the proper variables and reference axes be used for 
integration. To illustrate simply, Vy, should be integrated to obtain A and then 
Y'ce - ^ ; if Vy. is integrated to yield y* sg directly, then Vi becomes 

the difference of two large numbers and its accuracy may be reduced by several 
orders of magnitude. The above case is trivial, but the computation of the 
aircraft's translational velocity from the accelerations presents significant 
problems. If body axes are used, then during re-entry where the velocity ( V ) 
is large, rigid body rotations will cause large changes in the components of 
velocity which consequently may cause substantial errors in the actual accelera¬ 
tions (i.e. , a. x « Cl ). Howe* has suggested the use of the "H-frame"- X M 
directed along the projection of the inertial velocity on the plane normal to 
the earth radius vector, along this vector-for integration of the accelera¬ 

tions to obtain accurate velocities. The H-frame serves admirably for computing 
near-orbital velocities, but it is subject to difficulties at low speeds such as in 
landing. Here the inertial velocity is large compared to the aircraft's relative 
velocity ( V« r = 1500 fps as compared to V £ — 200 fps) and could cause errors. 

It appears that the best frame would be the T-frame- X r along V T and coin¬ 
cident with -so that the relative motion with respect to the earth would be 

computed with maximum accuracy. The use of either the H-frame or the T-frame 
requires extra transformation (as compared with body axes), but where compu¬ 
tational accuracy is limited (such as on an analog computer, and perhaps also for 
UDOFTT which has only 20 bit words as compared with the 32 bit words of the 
IBM 704) the use of such frame may be mandatory. No difficulty was encountered 
from this source in the present studies which employed an IBM 704 for solving 
the equations. For example, in the computation of the path of a highly eccentric 
orbit with rigid body rotations, the total error in position was found to be only 
30 ft in roughly a quarter of an orbit (30 ft out of 3, 000,000 ft). 

7.4 AERODYNAMIC REPRESENTATION 

Present knowledge of the aerodynamics of re-entry aircraft is not adequate 
to formulate an accurate description at this time. The best that can be done 
is to base the formulation on present estimation methods, but with the recognition 
that this formulation may be substantially in error and may neglect important 
variations which may not be uncovered until such aircraft are actually flight tested. 
However, it may be stated with confidence that for normal re-entry flight, an 
automatic control system will be needed and that direct manual control will pose 
a serious piloting problem and will only be used as an emergency'measure. 

Several possible simulations can be considered: simulation of the trajec¬ 
tory only, simulation of the attitude control problem only, and complete simu¬ 
lation. How each of these can be handled is described in Section 4.4.2, but 
the results can be summarized: trajectory simulation requires lift and drag 
data only, attitude simulation requires full description of the forces and moments 
but not necessarily as a function of flight condition, and the complete simulation 
requires including the dependence of forces and moments on the flight condition. 


* Fogarty, L>. E. and R. M. Howe: Analog Computer Solutions of the Orbital 
Flight Equations. University of Michigan, College of Engineering, sponsored 
by Air Force Contract No. AF33(616)-5864, to be published. 
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A recommended procedure for handling the aerodynamic data and arriving at 
suitable equations is presented in Section 4.4. The general equations for repre¬ 
senting the forces and moments in terms of aerodynamic coefficients are given 
by Equation 4. 12. As for the coefficients: basically, the longitudinal forces 
and moments ( C L , C D , ) require nonlinear representation as functions of 

attitude and flight condition (# , 0 ,S a >/V , , A ); the lateral forces 

and moments ( Cy , C 4 , C* ) can be represented by linear stability derivatives 
with nonlinear corrections added at hypersonic speeds. These concepts are 
summarized by Equations 4. 16, 4. 17, 4. 22 and 4. 23, and an illustration of 
the techniques required is presented in Figure 4.4. 

Normally, the force and moment coefficients are represented as function 
of & and # . However, the data discussed and equations presented in Section 

4.4.1 show that the direction cosines u.^/ V A , V A , u> a ./V* could be 

used to much better advantage in hypersonic flight. 

7.5 REQUIRED EQUATIONS 

In the preceding sections recommendations have been made concerning the 
geophysical and aerodynamic characteristics that should be included in the equa¬ 
tions of motion, and the coordinates and axis systems that should be used for the 
equations of motion. Here recommendations are made concerning the equations 
that are needed for the simulation. Since the degree of sophistication needed in 
the equations depends on the extent of the simulation involved and the task that 
is required of the pilot, three general types of re-entry simulation are consid¬ 
ered: trajectory and energy management simulation, short-period dynamics and 
handling qualities simulation, and complete simulation. 

7.5.1 Re-Entry Simulation 

Because trajectory and energy management simulation requires accurate 
geophysical models even when the navigational aspects are not considered, the 
least accurate representation that can be used must include a spherical rotating 
earth with an inverse-square gravitational force. If accurate navigation is 
required, then the oblate shape of the earth and its atmosphere, as described 
in Appendix C,must be included. The lift and drag ( C L and C D ) of the vehicle 
must be accurately represented, and also the side force (Cy ) unless it can be 
presumed zero. However, the attitude dynamics of the vehicle may be greatly 
simplified. In fact, if positive attitude control is provided the pilot, that is, 
if <k and 0 are made functions of control position only ( /& = 0 is presumed) , 
then the equations of motion can be reduced to the three degrees of freedom 
required to solve for the vehicle's trajectory. In this case the T-axis frame 
described in Section 7.3 would probably be the best frame to use for the equa¬ 
tions of motion. The body axes would only be used to define the attitude of the 
vehicle, for the pilot's attitude displays and for resolving lift and drag. The 
vehicle attitude angles can be resolved directly given Ot and 0 . 

Short period and handling qualities simulation requires accurate repre¬ 
sentation of the aerodynamic forces and moments with respect to all three axes 
to handle the lateral and longitudinal motions of the aircraft. In addition control, 
system simulation should be accurate. Since the emphasis here is on short 
period, the geophysical simplifications Jl - 0 and constant g as described in 
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Section VI may be conveniently used. Table 6.2 shows their negligible effect, if 
any, on the short period dynamics. Also, the mathematical simplifications, 
described in the next paragraph,may be used (though ?e, « 0 already if JT = 0). 

Complete re-entry simulation requires the complete equations derived in 
Section III except for the following modifications. The oblate shape for the earth 
and its atmosphere, described in Appendix C, must be included in the equations 
if accurate navigation is a requirement. The ( v er .) £ term in Equation 3.20a 
(follows Equation 3.72; term also found in Equation 6. 1) may be omitted. Further¬ 
more, basjd on the study_of the effects of this term, it can be concluded that 
(W)^ = (W) ss + tdgs-g x IV can also be omitted in Equation 3.20a. Therefore, 
the matrix equation 3.73 may be simplified to 

M ‘ M *] [ l S-Ss] {[V*r] + [*]} "W < 71 > 

thus affecting a considerable economy in the solution for & and /& . These 
conclusions assert the validity of the fourth assumption in Section 3.2.3 concerning 
relative air motion. However, these simplifications have been shown to be valid 
if [VJ , # . and /S are computed from Equations 3.68 through 3.71. If [/*]» 
i f , and /g are computed by integrating [V^ , OL , or yS directly, then these 
simplifications may not be valid. Since Equations 3.73 through 3.75 add consid¬ 
erable complexity to the over-all equations, their elimination is worth consid¬ 
ering. They could be dispensed with by numerically differentiating A and y& in 
the computer to obtain d and /& . Also, if & and yS are needed only for the 
aerodynamic representation, it may be acceptable in some cases to neglect their 
effect entirely. 

7.5.2 Landing Problem 

When the vehicle reaches its approximate landing location, it would be 
desirable to refer its position to some axis system fixed with respect to the 
landing field and to consider the earth flat. Longitude and latitude are no 
longer significant to the flight problem,and V? , 1$, , and V r can be integrated 
directly to obtain the position of the vehicle in rectangular coordinates. For 
simulating the approach and landing problem, the equations of motion can be 
greatly simplified by neglecting earth rotation, considering the earth flat, and 
assuming the gravitational force constant both in magnitude and direction. 

With these simplifications, the equations of motion revert back to the equations 
normally used for conventional airplanes. 

7.5.3 Computer Requirements 

It may prove of value to know how the IBM 704 storage capacity was 
allotted in the computer program used in this study. Of the total 3310 words 
used for the basic program, the breakdown is as follows: 

1. The input, including initial storage of all data, requires about 
710 words. 

2. The program itself requires about 2150 words of which standard 
subroutines comprise 1030 words. 

3. The output requires 450 words. 
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Included in the 2150 words of storage used for the program itself are an estimated 
550 words needed for the aerodynamic representation given by Equations 4. 16 
and 4.17. 


The time required for the time history solutions of Sections V andVI aver¬ 
aged about 40 computation cycles (time history points) per minute. 

7.6 CONCLUSIONS 

In conclusion, the results of this study show that significant simplification 
of the equations of motion can be allowed in their application to flight simulators, 
particularly in the calculation of the accelerations of the vehicle with respect, to 
the air. Additional simplifications can be made if accurate trajectories are not 
required, and certainly in some cases still other simplifications, such as small 
angle ones, can be made which have not even been considered in the present study. 


Because of the wide range of frequencies involved in simulating the full 
flight dynamics of piloted hypersonic or near-orbital vehicles, and because of the 
difficulties caused by such a frequency spread (Appendix D), it appears that new 
mathematical techniques for numerical integration, of adequate accuracy and 
favorable properties, need to be developed. If such new techniques are not 
forthcoming, it appears that considerably faster digital computers than are now 
presently available will be needed if the full flight equations are to be solved in 
real time. Several expedients can be adopted to solve the equations on present 
computers such as UDOFTT. The full flight of the vehicle can be subdivided, 
say into boost, orbital, re-entry, and landing phases, so that much simpler 
equations can be used in each phase than would be needed for the full flight. 

Also high frequencies associated with control systems, structural dynamics, 
etc. , can be neglected so that a long time internal can be used in the numerical 
integration. Analog computers do not offer a good answer to the problem because 
of the difficulty in obtaining adequate accuracy, and the very complex functions 
that must be generated to adequately represent the aerodynamics. Perhaps 
a solution can be arrived at by combining analog and digital techniques, using 
digital simulation where accuracy, complex function generation and logic are 
required, while using analog simulation for the high frequency elements of the 
problem. 
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APPENDIX A 

SUMMARY OF EQUATIONS OF MOTION IN ALGEBRAIC FORM 


This appendix presents the equations of motion as developed in 
Section 3.4 but, unlike that section, the equations are presented here in fully 
expanded algebraic form. Also, certain transformation matrices are 
included which are used .frequently in the equations. Each equation is given 
the same number which appears with its matrix form in Section 3.4. 

A. 1 EQUATIONS OF MOTION 

A. 1. 1 Equations of Motion 


U « -2- + - fur 

-w* 

sis' • ~j-- - ru + •jDur’ (A. 37) 

^ * -J- ++* 


■ j + . /** , "t * lihfj [( z » rt\ 

J * - J *t I * 1 } J *; 

+ ^7 [fe - x r) ** '*] 

f '- 1 * 6*'V*)] (A. 38) 

* * j f' *, j * M * *1 ~T*-I 1 [A “* 1 #^ ^ J *»/ 

x * ' x *t x * J f 

* - /*, . [(*« - v** -■*■«> H 

*n X J *XJ 


A. 1.2 Trajectory Kinematics 

(«#• a ess r)+ s* (*•* 0 •*» -«'*» * tM 

+ ur(c»9 9" c#* 0 •<> Q *• *1* 0) 

V x = u(c 0 »e si»r)+'*'(*• t>n 0 +0) (A. 41) 

■f ttrfsJm f o** 0 •/" & ~ ‘’O* f sin 0) ~ ^e»» 

« af- sin ur(cm% 0 co*o) 


t 
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v *r 3 • /l >V # «« 

(A. 42) 

l/ r - {v* + Vf V T >0 

(A. 43) 

% ■ fi? *iv"i li'O 

(A. 44) 

* j v A * * - - j 

| (A. 45) 

. os ^>o 

' 02^ >-f*V#r* 

(A. 46) 

< s * 

(A. 48) 
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-fs> t - jo -SI (to* ■*> cos "it coo 0 *■ sin i> sin 0) 

<f. g » +SL (to* ■$> sin It cos 0 - a os H cos ft sin & sin 0 *■ sin cos 9 sin 0) 
)\ ■ r -SI (to*f> sin 0" sin 0 y- cost co* 0- sin & cos 0 - tin i> cos 9 to* 0) 

(A. 57) 

Details of the vehicle orientation angle switching mechanism are 
given in Section 3.4.3. 


^S • \ COS 6i I r * \ COSG^ 
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c**e t 


+ i sin #4 


(A. 55 


A.1.4 Relative Wind Kinematics 


U^m u - Wf, (cos Geo* ft) ~ (V mr + Wx )c*s 6 sin ft • tin 0 - M x 

sir^ * sv - Wj (sin 0tin 0 co* ft - tin ft to* 0) 
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/& - rj— (~ sin /4 <»» » + c»t/9 - < 4 ^ M/n/t o/n m) (A. 74) 
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i +JL (cos P sin If cot 0- too pcot ft sin 9 sin 0+ tin P cos Q t in 0) 
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(A. 76) 

A. 1.5 Forres and Moments 


■x/’ 1 '* se « 

5£ » (A. 79) 

^'X/O 3 «a 

-i^y'skcj 

.+/> IJ* S«C„ (A. 80) 

F h * ***• {* 5 “") t,m 0 C0 * e (A. 81) 

% * ’"*(£) 

Equations A. 37, A. 38, A. 48, A. 49, A. 50, and A. 55 or A. 55s are simul¬ 
taneously integrated using numerical methods. 
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A. 2 TRANSFORMATION MATRICES 
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APPENDIX B 

EXTENSION OF THE EQUATIONS TO LARGE 
THRUST AND VARIABLE MASS 

Equations of motion derived in Section III of this report are, restricted 
to a rigid, constant-mass, re-entry glider. In this appendix equations are 
developed to describe the motion of a booster-glider combination. The general 
vector equations of motion for a jet-propelled body of varying mass are pre¬ 
sented in integral form in Section B. 1. Application of the general equations 
is outlined in Section B.2, where the various integrals are evaluated by making 
use of reasonable approximations, and a variety of problems related to use 
of the equations are discussed. 

B. 1 THE.GENERAL VECTOR EQUATIONS 

In Reference B. 1* there are presented the full details of a derivation 
of the general equations of motion for a jet-propelled body of varying mass. 
Here in this appendix only pertinent assumptions, final results, and a mini¬ 
mum of detail are presented. 

A body of constant volume T whose surface area S has fixed size 
and shape is assumed**. Mass density /O is variable within the body, a 
function of time and location. Associated with the body is a reference frame 
X , y , £ which is fixed relative to the body and with origin at some arbitrary 
point B . This frame will later be identified as the body axis frame of 
Section III. 


The linear and angular momentum equations of fluid mechanics can be 
applied by considering the surface of the body to be a control surface encom¬ 
passing a region of space, the body itself: 


E ? + u + 00 * X + *Ai)(' 1/ k •»)*<* (B. 1) 


E^r /° [^/s *(V+£> k 

+ (B.2) 

With the exception of definitions given below, the notation is consistent with 
that elsewhere in the report: 

£ position vector of any particle relative to body axes, 

» (i) a velocity vector of any particle relative to body axes, 
5f unit normal vector to surface S , directed outward, 

* Other sources of similar, but not so general, derivations are listed in 
Reference B. 1. 

** Although these assumptions are not strictly essential, they afford con¬ 
siderable simplification in the derivation. 


122 



9 


NAVTRADEVCEN 801A 


£ F summation of all external surface and volume forces 
i acting on the body, 

£ summation of all external moments on the body about 
the origin of the inertial frame, 
differentiation relative to the Inertial frame, 

- (•) time differentiation relative to body axes, 

* J* ( volume integral over the volume contained by S , 
j Bur f ace integral over closed surface 3 . 

In the terminology of fluid mechanics, the volume integral in each equation 
is called the "nonstationary" term and the surface integral is called the 
"convective" term. 


By using the definition of the mass of the body 
7K 9 

and the definition of net mass efflux 


lux 

0 - f /°("*-S)d* 


(B. 3) 


(B. 4) 


Equation B. 1 can be simplified somewhat to read: 

i -~T + fc). * + jf fa * *)** 

f *i fa -n)ttd ^ V- 2 *> Jtj f> r.(*r M n)d+ 


Y.P ^ 

dt 


(B. 5) 


In this equation, the following new notation is used: 


position vector of center of mass relative to body axes 

position vector of any particle in the body relative to the 
center of mass 


These vectors and some of those defined previously are depicted in the sketch 
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The summation of forces on the left-hand side of Equation B. 1 includes 
many surface and volume (body) forces of which only gravity and surface pres¬ 
sures are considered here. The forces can be written as 

JcLT + ^ =£ T (B-6) 

in which these definitions apply: 

^ gravitational force per unit mass (acceleration of gravity) 
which is assumed to be uniform throughout the body 

p pressure acting on outside of surface S . 

The pressure, -p , can be broken down into two parts: 

P - Pa * (i° ~ Pa) 4 -Pa**** (B.7) 

in which p is the ambient (atmospheric) pressure assumed uniform every¬ 
where in the immediate vicinity of the body. Thus, since is uniform over 

S , the surface integral over S of ^ is zero, and Equation B.6 
becomes 

•h*J v- I - ApJt d<? (B. 8) 

In Equation B.8 it is possible to identify the integral of bp over S as the 
aerodynamic force on the body provided other surface forces, which are implicit 
in the integral, are zero. Further evaluation of the pressure integral requires 
specific information regarding the vehicle to which the integral is applied. 

The result of the evaluation of Equation B. 1 is 

' ^7?] * **(%)**} 

*[£i 

Equation B. 2, the moment equation, can be made more meaningful 
by modifying it so that the summation of moments is about the body axes 
origin, £ M 0 ’ 

Z"* -Z^t * Z^ . 

ujLj |> g (v+£3* dX f V*j /o(7+ S>xK+~j,)dr 

*'j /> [*»{?*■*»* ’*)*<* (B. 10) 

Following techniques analogous to those used with the force equation and using 
Definitions B. 3 and B.4, Equation B. 10 can be simplified: 


* (B. 9) 
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Eff. -15 *JL ,r 

* T [^*(*** *)](*Sk •*)<*<* + { fi (** '**)("* ■ *) d<f 

Js (B. 11) 

Except for the fact that /O is not constant, the first two terms are those 
which apply to a rigid body whose moment equation is evaluated relative to an 
arbitrary point (Reference 3. 1). The last three integrals, due entirely to 
internal mass motion and mass-efflux, would vanish for a rigid body. 

The left-hand side of Equation B.2 consists of the summation of moments 
due to all surface and volume forces acting on the body. As before, only the 
forces due to gravity and surface pressures are considered; in addition, it 
should be noted that ambient pressure produces no moment: 

F t x x J + r\ x £ - A f> * d<f + J x Ay? ft d& A (B. 12) 

Then, after combining Equations B. 11 and B. 12, the resulting moment equation 
is 


X 


*)](*„ ■ Z)d« 

(B. 13) 


The first term on the left-hand side of the equation is the gravitational moment 
about point B. The two surface integrals on the left-hand side would be the 
moment about B due to aerodynamic pressure if there were no other pressures 
(distributed or concentrated forces) acting on the surface of the body. 

B.2 APPLICATION OF EQUATIONS TO ROCKET VEHICLES 


In this section, a short discussion is presented to facilitate applica¬ 
tion of the equations of motion developed in the preceding section. A rocket 
vehicle is chosen as the vehicle of interest. It consists of a glider attached 
to two booster stages which eventually separate away from the glider during 
the initial phase of flight. 

The development of equations throughout this report has been directed 
toward a description of the glider's motion; the objective is the same here. 
Equations B.9 and B. 13 are written for the motion of the reference frame whose 
origin is B and which is fixed relative to a body of varying mass. If this 
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reference frame is identified as the body axes frame of Section III, the equa¬ 
tions of motion developed in this appendix describe the motion of the glider in 
terms of translation of its center of mass (B) and rotation about its center of 
mass. In addition to interest in the glider's motion, it is desirable to use a 
reference frame which is fixed relative to a body because forces and moments 
are more easily described. Specifically, the origin of the body axes frame is 
the glider's center of mass (assumed fixed); the X axis is directed toward the 
nose of the vehicle, the *1 axis is directed downward, and the tj axis is 
directed toward the right. 

B.2. 1 Evaluation of the Remaining Surface and Volume Integrals 

Equations B.9 and B. 13 contain surface and volume integrals which 
must be evaluated; otherwise.the equations are not very useful as equations of 
motion. Evaluation is made for a rocket-propelled vehicle. 

The assumed planar and circular cross-sectional area ( A % ) of the 
exit end of a rocket motor is the only area which need be considered in the 
surface integrals. No mass enters or leaves the vehicle elsewhere. It is 
assumed that svjg over the exit plane can be replaced with an effective value 
<*r e which is uniform in magnitude and direction. The vector 7* over the 
exit plane can be replaced with a_vector directed to the effective "center" 

of the exit plane and a vector Ar directed from the effective "center" to 
any other location in the exit plane: 


r* e *■ Ar 


(B. 14) 


7 l stn FLAM* 

The effective "center" is defined such that 

f p /* (sv-j, • x)dd - r c Q * [ p dr (^ M •’x)c(& » Q (B. 15) 

J A a \ 

However, when appears twice in a triple vector product (see Equation B. 13, 

first surface integral on right-hand side), the integral containing this product 
does not vanish. Reference B. 1 shows details of the evaluation for this case 
in terms of the radius of gyration of the exit plane. 

Equation B. 13 contains the volume integral, $ T pH which 

can be separated into two parts using 7^ and Fj , but this process leads to 
many small terms involving and its derivatives in addition to the integral 

dtSf./ 0 ^* *** dr which is really just as bad as the original integral involv¬ 
ing X . The original integral cannot be considered always unimportant because 
it includes the effects of internal rotating parts, fuel sloshing, etc. Without 
any knowledge of the internal motion in the rocket vehicle, nothing can be said 
about the value of the integral; therefore, it will be carried along as it is. 

The exit-plane portion of each pressure integral in Equations B.9 and 
B. 13 is usually considered part of the total thrust force or moment due to the 
rocket motor. It is assumed that the exit pressure (-p % ) of the rocket motor 
is effectively uniform over the exit plane; thus 


Jam 


(B. 16) 
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The remainder of each pressure integral is now defined as an aerodynamic 
force or moment, as the case may be. Evaluation is usually made in terms 
of experimentally or theoretically derived aeroydnamic coefficients (e.g. , 
see Section 4.3). 


Making use of the assumptions described in the three paragraphs above, 
Equations B.9 and B. 13 can be written in more meaningful forms. The |orce 
equation is: 

P 9 * ?+ + *T -X [(ft, * * * V ] 

+ >1 j/Hi * * **C +2* *(K)b +8*(&*K)]+ *(%) 8 

(B.17) 


in which 

K 



6pn d 


F t = - Q [/tr e +2 to * F. y ( \) a y r e Q + (-P A A m 
The moment equation becomes: 


(B. 18) 
(B. 19) 


B.Z.I'a.f'S (b .2 «) 

in which = "M. )\ xg 

* F c x I -Apfidcr + I r»- x(B. 21) 

Js-At 

M T = - Q x (to * n B t) + j- (u> + *>x * '*’•] 

* F es * (*?» “■*) *• (B. 22) 

Certain symbols and notation used in these equations must be defined: 

t* m position vector from center of mass to effective center of exit 
area 

unit normal to exit plane; positive outward 

effective value of over exit plane, or exit velocity; 

positive outward 

y* a9 position vector of effective center of exit plane relative 
to point B ( Y* at - F* + r* e ) 
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^moment of inertia tensor of the entire vehicle (booster and 
glider) relative to point B 

radius of gyration of exit area ( 4 e = radius of exit area) 

component of to in X body axis direction; e% is unit vector 
and Uf x> =-p , the rolling velocity 

I ( )d& surface integral over entire surface of vehicle except the 
exit area Ae 

B.2.2 Obtaining the Desired Accelerations 

At this point, the force equations would be solved for (^)g and the 
moment equation for & as is done with Equations 3.35 and 3.36; however, 
inspection shows that both (VJ a and S> appear in both Equations B. 17 and 
B.20.^ Therefore, Equations B. 17 and B.20 must be solved simultaneously 
for (V) B and S> . - 1 

_ Simultaneous solution of Equations B. 17.and B. 20 results in both F 

and A7 appearing in each of the equations for fl^and &> , which precludes 
the use of aerodynamic derivatives such as C L ^ , or any other 

derivative whose application depends somehow on the linear and angular 
accelerations ( and to ) . If the affected derivatives are modified to 

be used with the linear and angular accelerations explicitly (e.g., components 
, etc. instead of dt , , etc. respectively), then again Equations 

B. 17 and B.20 can be solved for ($) g and eo . This is, however, easier said 
than done: Vectors (VJ g and to are rates of change of inertial velocities rela¬ 
tive to the body axis frame whereas the pertinent aerodynamic derivatives are 
with respect to the rates of change of the vehicle's linear and angular velocities 
relative to the air, vectors and • 11 may prove more reasonable, 

if not necessary, to neglect all aerodynamic derivatives such as those mentioned 
above. In Section III of this report, it is shown that without special treatment 
such as an iteration, the equation for [v] (Equation 3.37) cannot contain 
derivatives with respect to d? and /& but, because of dt and /& do not 
depend on [&] , derivatives involving dt and /? can be used in Equation 3.38 for 
[a*] . The situation is worse in. Equations B. 17 and B.20. These equations 
must be solved for [v] and [o>] and, since both [F\ and [Afl appear in the equa¬ 
tion for [v] , even the moment derivatives with respect to dt and 
(which are more important than the analogous force derivatives) cannot be used 
without special treatment. 

B.2.3 Data for the Combined Booster-Glider Vehicle 


In Equations B. 17 and B.20, center of mass location ( r^) and the 
moment of inertia tensor (f') are required in terms of instantaneous mass: 


y'c • r c C-m) 

I’- T'(-m) 


(B. 23) 


The specific form of Equations B. 23 may be either functional or tabular and, 
in fact, reasonably smooth tabular data can be curve fit if desired. At any 
instant, mass is given by numerically integrating ~m - - : 
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* »€ 0 - I Q dt (B. 24) 

Jo _ _ 

Equations B.17 and B.20 also contain first derivatives of f| and I* which can 
be obtained directly by numerical differentiation with respect to time, or 
indirectly by either functional or numerical differentiation with respect to mass 
(depending upon the form of Equations B.23) according to the scheme: 


(K) s 

n 



(B. 25) 


The second derivative of center of mass location, (fy) L, can be easily found 
by numerically differentiating whereas a scheme similar to that shown 

in Equations B.25 is more difficult because * unlike F* e , is a function 

of both »e and Jn . Fortunately, the velocity and acceleration of the center 
of mass relative to the body are small; therefore, inaccuracy in the numerical 
differentiation is of little consequence. 


The quantities pertaining to rocket motor operation ( <5 = -"tit , Q , 

Ar-q , , etc.) can be related using equations governing the injection, com¬ 

bustion, and expansion of the propellant gases. However, a useful expedient 
for situations not requiring utmost accuracy is to treat the mass flow rate as 
a function of throttle position. 


Q * Q (* T ) 

Such dependency is useful if the rocket thrust magnitude is controlled by a 
pilot or, lacking a pilot, S r can be made a function of time or some other 
variable. The quantity & may then be obtained by numerical differentiation 
of with respect to time. 


Equations B. 18 and B. 21 represent integral forms of the aerodynamic 
force and moment respectively which, as stated previously, may be evaluated 
in any desirable manner, but most probably by using the method described in 
Section 4.3. It must be remembered that during the boost phase of the flight 
each separation of a stage changes radically the over-all aerodynamic repre¬ 
sentation of remaining portions of the vehicle; therefore, effects of separation 
must be considered in the evaluation of Equations B. 18 and B.21. 

B.3 REFERENCES 
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APPENDIX C 

EXTENSION OF EQUATIONS FOR EARTH OBLATENESS 
AND SURFACE ABERRATIONS 

C. 1 EXTENSION FOR FIRST-ORDER OBLATENESS 

It is assumed in Section III that a spherical representation of the 
earth is entirely adequate for purposes of re-entry flight simulation. This 
implies that both the geometric and gravitational effects of oblateness are 
negligible. Consequently, the simulation equations for earth size, shape, 
and gravitational field are fewer and simpler providing a substantial saving 
of digital computation time. This saving is certainly a desirable feature of 
the spherical earth representation. 

The models which constitute adequate representation of the earth's 
size, shape, and gravitational field are determined by the scope and purpose 
of the flight simulation. Simulation of the boost and orbital phases of a 
vehicle such as Dyna-Soar probably requires an oblate representation of the 
earth, while the re-entry phase can be adequately represented using a 
spherical-earth model. In the event that the more precise representation 
may be needed, the process for inclusion of earth oblateness with the equa¬ 
tions of motion is developed here to first-order accuracy. 

Inclusion of earth oblateness should generally result in both geomet¬ 
ric and gravitational effects although it is not really necessary to include 
both. The development in Section C. 1. 1 and C. 1.2 implies that the two 
effects can be separated for simulation purposes if desired. The geometric 
effects are, roughly speaking, larger than the gravitational effects so that it 
might be desirable to include the former but not the latter. 

Geometric effects can be arbitrarily divided into three groups: 

1. Location of vehicle-latitude, longitude and altitude 

2. Atmospheric models-density, temperature, etc. 

3. Vehicle orientation-the Euler angles 

Longitude is unaffected because the oblate shape of the earth is axially 
symmetric. Latitude is different by no more than about 0.2*,. a small angle; 
but in terms of distance on the earth's surface.this corresponds to about 12 
nautical miles. Altitude is much different from that of the mean spherical 
earth. For a given constant distance from the center of the earth the oblate 
earth causes the altitude to be greater by about 47,000 ft at the poles and 
smaller by about 26,000 ft atthe equator when the spherical-earth radius of 
Equation 2.3 is used for comparison. If altitude above the oblate earth is 
used with the models of atmospheric properties given in Section II, these 
properties become "oblate" to conform with the earth's shape. The Euler 
angles are affected only slightly, the largest effect being less than 0.2* for 
'any of the angles; but, since these angles are used in coordinate transforma¬ 
tions, the slight effect should be considered important. There are numer¬ 
ous sources for further information on the effects of oblateness; References 
C.4, C.5, and C. 6 are examples. 
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C. 1, 1 Earth Size and Shape 

The earth is assumed to be oblate spheroidal; that is, it is an ellipsoid 
of revolution whose radius ( ) at any geocentric latitude ( ^ ) is (an alter¬ 

nate form of Equation 2.2) 

a(/~e) 

" //- (Z6-6*) co,**' (C * l) 


For this equation: (data from Reference 2.3) 

CL , the equatorial radius = 6, 378, 163 meters = 20, 925, 732 ft 
£ , the flattening = ^ 

A first-order approximation in £ for can be obtained in the form: 


X00 « *('-£ 2 *0 


(C. 2) 



The quantities to be obtained are -ffc and line CB - -A> r in terms of , 

CL , and € which are known. This is done by developing corrections 
for and »C , the variables ordinarily computed. Reference C.2 supplies 
expressions for - P * £ and < T but using notation different from this 
appendix. Equation 4 of Reference C.2 is a series expansion for f which 
is modified and simplified to result in the equation: 




S/> t t> 


(C.3) 
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This expression has an error of less than ±0. 5% when compared with the 
more precise expression in Reference C.4 and should, therefore, be adequate 
for simulator use. 

Refer to the sketch above. Line 1T& is the true altitude, "A r . The 
radius of the oblate earth at geocentric latitude, *> , is 35 = y? ## . Altitude, 
&8 - A , is computed from the equations of motion of Section III. The 
altitude along the position vector 08 is 53 = JTS - JZs . It can be 

shown, by comparison with the altitude computed from Equation 6 of Refer¬ 
ence C.2, that 5"3 ~ A r with an error less than ±0. 1%. Now, this means 
that a correction, £3 = 4 £ , which is the difference between the oblate- 
earth radius and the mean spherical-earth radius (as used in Section III) can 
be added to the altitude, A , to obtain the altitude, A r , above the local 
surface of the oblate earth. The following equations define this process: 

&< = *oo (C.4) 

AA-(A-K 0 )-*,6si»*f> (C. 5) 

A r c A + &A (C.6) 

The mean spherical-earth radius is given by Equation 2.3. 

Methods for correcting spherical-earth latitude and altitude for oblate¬ 
ness are discussed in Section C. 1.3. Further information on the problem of 
vehicle position relative to the oblate earth is presented by Reference C. 1 
in a development similar to that given here. 

C. 1.2 Gravitational Acceleration and Force 


The gravitational potential function is given by Equation 2.5, but is 
repeated here for convenience: 


u = 




(C.7) 


in which P^ is the n 1 * 1 degree Legendre polynomial in sin . 

a = 6,378, 163 meters = 20,925,732 ft 
= 1 
j; = o 

J z = -1082.61 x 10' 6 
J, = +2.05 x 10' 6 
= +1.43 x 10' 6 

Gm m = 3.98602 x 10 14 meter ^ = 1.40765 x 10 16 ft 3 /sec 2 
• sec* 
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The second degree Legendre polynomial will be used later; it is: 

P * *7 (3 *'"**>-/) (C.8) 

A listing of Legendre polynomials and their properties is given by Jahnke 
and Emde, Reference C. 3. 

The gradient of the scalar potential function is the vector accelera¬ 
tion of gravity, thence the force is given by the product of mass and accelera- 


tion: 

9 • V V 

(C.9) 


r s • ** V V 

(C. 10) 


The operator (7 must be expressed in terms of a particular coordinate 
frame, here the spherical# coordinate frame is chosen: 


r 7 a L 3 tj- _/_ d — 9 _. 

9 ^ ** * 9X ** " (c U) 

in which , 8^ , and !,» are unit vectors in that reference frame. 

The components of are then obtained in the same frame. 

Using only the second harmonic which is adequate for simulation, the 
following matrix equation for components of JL in spherical axes results 
from Equation C. 10: 9 



Here the matrix [VU] represents the components of the vector vu in the 
spherical coordinate frame. Equation C. 12 should be used to replace the 
corresponding equation of Section III. 


# In this appendix, the term "spherical" denotes the frame previously called 
"earth-surface" in Section III. This usage is explained in the fourth para¬ 
graph of Section C. 1.3 which follows . 
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C. 1.3 Modification of Equations of Motion to Include Oblateness 

This section discusses ways to incorporate earth oblateness correc¬ 
tions to the equations of motion developed in Section III. Certain concepts 
presented in Sections C.l. 1 and C. 1.2 are used in addition to, or to modify, 
the equations of Section III. The result is an adequate first-order representa¬ 
tion of oblateness. 


What variables are used? The vehicle's six inertial velocities ( U , 
tUr , £ , y ) are retained. The variables A , X are essentially 
spherical position coordinates which, of necessity, are retained. Relative 
to the local surface of the oblate earth, the vehicle's angular position is 
given by the new Euler angles (0g , ). These twelve are the variables 

chosen for simultaneous numerical integration. 

Position of the vehicle relative to the oblate earth is described by the 
pseudo-spherical coordinates (A , X , P ) appropriately modified. True 
altitude (A r ) is obtained by adding the correction (4 j€) of Equation C.5 to 
the given altitude (A ). In like manner, geodetic Latitude (^ ) is obtained 
by adding the correction ( 5 ) calculated using Equation C.3 to the geocentric 
latitude (P ). Longitude ( X ) is unaffected by oblateness. Thus, vehicle 
position is simply handled. 

The earth-surface ( ES ) reference frame has been previously defined 
as the frame whose axis is directed downward and normal to local 

(smoothed) surface of the earth; axis is directed northward; axis 

is directed eastward. The plane is then parallel to the local earth's 

surface. In Section III the earth-surface axes were associated with a spheri¬ 
cal model of the earth,while in this appendix the earth-surface axes are 
associated with an oblate-earth model. Consequently, the earth-surface 
frame here is not the same as that of Section III. Since certain variables in 
this appendix must be expressed in spherical coordinates, the earth-surface 
(ES ) frame of Section III is retained, but denoted differently as the spherical 
( SPH) frame. Then, the frame associated with the oblate earth is denoted 
as the earth-surface ( ES) frame to be consistent with the definition of that 
frame given above. The origin of each frame is the same point (B), the 
center of mass of the vehicle, and the axas collinear with the 

axis. This new notation should be carefully observed throughout the remain¬ 
der of this appendix. 


Transformation from the spherical (SPH) axes to the earth- 
surface (£■$) axes is obtained by using 


M- 


/ 

o 


-« 


* K 

f o 

0 / 


(C. 13) 


Here the small angle approximation for sines and cosines is valid because 
f is always less than 0.2*. The complete transformation from earth (E ) 
axes to body ( 3 ) axes is 
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[ L *-t] ■ [*-,-) [w] - H [«.] [«] [«] M [x'jM < c - 14 > 
where * M [«.] [«■] [?] ( c - i5 > 

The transformation here defined is identically the trans¬ 

formation of Section Ill except that here it is denoted differently and consists 
of four consecutive rotations. A further definition is made for convenience. 

[*«-«]■[*«] [».][«] « c - i6 > 

This new matrix transforms vector components from the earth- 

surface frame to the body frame. 

The velocity of the vehicle relative to the earth must be calculated in 
two ways. First, in order that the pseudo-spherical coordinates may be 
obtained, the vector [Vj is calculated in spherical ( SPH) axes. The .com¬ 
ponents of this vector, UV g ) t0 /\ . are required for computation of 4 , 

X , and from Equations 3.48, 3.49 and 3.50 respectively. The magni¬ 
tude ( Y % tB ) of the position vector, needed for Equations 3.49 a*hd 3.50, is 
computed by adding altitude ( A ) to the mean radius (^) of the spherical 
earth. Second, the relative velocity in earth-surface axes is 

obtained as 

[flfci,] * [r] [ftbJu] <c. i 7 ) 


The Euler angles (££ , 9 0 , 0 0 ) are obtained by integrating an 
equation identical in form with Equation 3. 55: 



\*\ 


' o' 


* 

■M M - [•%] 

■f 


in which 


M 


[*.} 


i/co s 9 0 O 

O / 

S/'/f /C0» O 0 O 

‘//cos O 

O / 

sin 9 0 /cot o 


f] N>.]' 

o] [«] [«] MM 
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The rate of change of geodetic latitude ( ^ ) is obtained by directly adding two 
collinear angular velocities: 

H * t (C. 19) 

for which differentiation of Equation C.3 gives Jf : 

2 P 64 * 2 i> _ £ s/n 2 P 
' ft,* 



(C. 20) 


The vector is given by Equation 3.57 except that [^-0-f] is now defined 

as in Equation C. 14. 

Velocity of the vehicle relative to the air is calculated in essentially 
the same manner as in Section III: 

M = M -K-«] {M * M}-W (c - 2l) 

Here, the wind velocity vector [iv] has components in the new earth- 
surface frame which are calculated using geodetic position variables such as 
and ^ . It should be noted that [y.r] the same in either earth-surface 
or spherical frames. If the gust velocity vector [/•<] should depend on the 
vehicle's position, it must also be calculated using geodetic position coor¬ 
dinates. The components of [V*] are used to calculate its magnitude , 

and dC and /ff 

To obtain [V*.], Equation 3.73 is written in terms of the oblate- 
earth quantities already defined in this appendix: 

[i ] a M * h* *] [*--•] {M * M} - [>] 

- [(**”-*)**,*][&] M'W [^ er i] +[%\ M} 

(C. 22) 

Then i and /9 are calculated by using Equations 3.74 and 3.75 with the 
components of [i*,] and other quantities associated with the oblate earth. 

In addition to the wind and gust velocities, other geophysical quantities 
are calculated as functions of geodetic position variables. Specifically, atmos' 
pheric properties, such as density, temperature, and pressure, are functions 
of A r , the true altitude above the oblate earth. 

Gravitational force for an oblate earth is given in terms to the second 
zonal harmonic by Equation C.7. This force, in terms of spherical ( 5PH ) 
coordinates, must be transformed to the body axes coordinates: 

[*] * M ["sU ,c - 23 ' 
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where again is defined by Equation C. 15. Aerodynamic forces and 

moments are calculated using oblate earth quantities V x and fi (A r ) • 

C. 2 SPECIAL SUF FACE ABERRATIONS 

Certain small-scale geometric aberrations of the earth's surface such 
as mountains or mountain ranges, etc. which would comprise obstacles to 
low altitude flight, near landing for example, can be simulated. One possible 
method for providing such obstacles is given here and another more precise 
method is mentioned. 

The equation 

= «..• j--: -ii)utjr x (x~-x) (c.24) 

with the conditions 

I " *%\ < 

I ^ I < 

defines what can be called a "cosine ridge" because it can be longer than wide 
or vice versa depending on choice of Lp and L x . Here: 

is the west-to-east dimension of the ridge in 
radians on the earth's surface, 

L x is the south-to-north dimension of the ridge in 
radians on the earth's surface, 

H is the height of the ridge in feet, 

•/> , A__ are latitude and longitude of pole or center of 
** ridge, 

is the correction to the vehicle height in feet. 

There are no restrictions on the magnitudes of Lp , L x , or H . Also, 
geocentric latitude () can be used in place of geodetic latitude () if 
the earth is assumed spherical as in Section III. 

If a number of separate ridges are placed near one another, they 
combine to produce a "mountain range". At any particular location, the 
total height correction for the vehicle is obtained by summing the cosine 
ridges, including the conditions on each A as shown above. The true 

height of the vehicle above this obstacle is given, then, by the addition of 
the summation to Equation C.6. 


A\ r ■ A + &A 




(C. 25) 
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It must be remembered that the true altitude to be used with the air 
density, air temperature, wind function, etc. is that given by Equation C .6 
without the correction due to the obstacle. The altitude, , as obtained 

from Equation C.25 should be used only as a cockpit display or as a program 
stop if it should become negative, indicating a crash. 

The computation of A from Equation C.24 can be facilitated by 
using the first two terms of the power series expansion for the cosines. 

This is allowable because the obstacles need not have any particular repre¬ 
sentation. In fact, the cosine representation was used only because it 
resulted in an obstacle with roughly the appearance of a "mountain" or 
"ridge". 

There are methods for obtaining a more precise simulation of a 
portion of the earth's surface, if such is desired. For example, a two- 
dimensional array representing elevation ( £ ) of terrain at intersection 
points of grid lines of longitude ( A ) and geodetic latitude (#£ ), or 
equivalent coordinates on the earth's surface, can be stored in the memory 
of a simulator such as UDOFTT. For a given position of the vehicle, the 
computer interpolates the data of the array to obtain elevation, £ ( A , P* ), 
which modifies the vehicle's altitude. Using an oblate earth model, within 
the region covered by the array, the true altitude of the vehicle is: 

A t = <+fit.- (C. 26) 

It is possible to use orthogonal functions in two independent variables to curve- 
fit elevation contour data and then use the orthogonal functions in the simulator 
computer. However, this method is probably more costly in storage and 
computation time than the interpolation method above for equivalent accuracy 
of simulation. 
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APPENDIX D 

SUMMARY OF DIGITAL TECHNIQUES USED IN PROGRAM 

In the sections which follow, the important digital techniques used during 
programming and solution of the equations of motion are discussed. Needless 
repetition, is eliminated, however, so that only topics not covered elsewhere are 
treated in detail here. The three main section are devoted respectively to 
general mathematical techniques, numerical integration, and obtaining solutions 
of the equations. A simplified diagram of the computer program (Figure D. 1) 
is included to illustrate the computation scheme. 

D. 1 GENERAL MATHEMATICAL TECHNIQUES 

Many of the digital techniques mentioned in this sectic l are discussed in 
Reference D. 1 which is a general description of the IBM 704 program evolved in 
this study. The reference contains a description of the characteristics of the 
problem, a mathematical statement of the problem, the method of solution, a 
series of flow diagrams of the computer program, and instructions in the use 
of the program. Because they are treated in Reference D. 1, the computer 
routines are discussed only briefly in Section D. 1. 1. 

The remaining sections, D. 1.2 and D. 1.3, deal with axis switching and 
the control system respectively, both of which are described in detail elsewhere 
in the report; therefore, the two sections will be short. 

D. 1. 1 Routines for the Computer 

Of the digital routines which facilitate use of the program, the input and 
output routines are among the most important. The input routine provides for 
loading of any storage location in the program, consequently it facilitates 
changing of conditions for different solutions. Input and output routines both 
use floating-point decimal notation. Specifically, the output routine prints 
floating-point numbers which at any given time may be used directly as initial 
condition inputs for a continuation of the same solution or a different solution. 

A time history solution may be halted upon reaching some condition of a com¬ 
puted parameter, or upon reaching a preselected time. To increase the computing 
speed, 900 core storage locations for each output variable are used to store 
results so that mechanical operations (of tape feed mechanism) .are not performed 
after each computation cycle. The print-out from core storage to tape is made 
at one time when the entire 900 locations become filled or sooner if desired. An 
option exists for storing only the solution of each n^ 4 computation cycle. The 
input routine code number is RS 0046, and the output is GL OUT 2. These are 
SHARE sub-routines. 

Program outputs require three magnetic tapes, each handling about one 
third of the total number of variables. Arranged according to the tapes in which 
they appear, the variables are: 

1 . U , ^ , ur ; ^ , r ; 

2 . 0 , <9 * , M t 

3. (t , /S , ? , V* , IV , V e , r r 
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where Cfj = -j p V* and M r is the moment about the tf body axis due to 
the controller described in Section 5. 1. 6. 

Other standard SHARE sub-routines are used for square root, sine, 
cosine, exponential and natural logarithm. These are listed with their code 
numbers in Reference D. 1. Because of errors arising from truncation of 
numbers by the computer, a SHARE sub-routine for arcsine (CL ASC 1) was 
modified so that the program would error stop if the argument of the arcsine 
were greater than 1.001, or the argument would be used as unity if it actually 
were less than 1.001. 

Three matrix sub-routines were developed ' specially for this IBM 704 
program. They are routines which: 

1. premultiply a 3 x 3 matrix by another 3x3 matrix, 

2. premultiply a 3 x 1 (column) matrix by a 3 x 3 matrix, and 

3. obtain the transpose of a 3 x 3 square matrix. 

A diagnostic routine called FID© was used extensively during checkout of the 
program. This routine enables one to obtain a print-out of the content of any 
storage location used during the computation cycle. 

D. 1.2 Axis Switching 

Vehicle orientation may be approached from two points of view: the use 
of direction cosines and the use of Euler angles. The latter have been chosen 
here because they provide a more meaningful physical picture of orientation. 
The Equations A. 55 for £ and P contain a (cos Q ) _1 factor which becomes 
zero as & approaches *90*. An actual re-entry vehicle can possibly operate 
at this condition; therefore, the program would be severely limited if it could 
not adequately handle the Q —* *90* case. This is the purpose for the axis 
switching mechanism described in Sections 3.3.2 and 3.4.3 and in Reference 
D. 1. 


A present restriction on use of the program is that the vehicle not be 
allowed to pass directly over either pole of the earth. This is a result of using 
latitude as a position variable for the vehicle. Note that Equation A. 50 is not 

defined when cos P approaches zero; that is, when P -► *90*. A switching 

mechanism similar to that used for the Euler angles could be developed to 
overcome the polar flight restriction. 

D. 1.3 Control System 

The basic digital technique use of the controller described in Section 
5.1.6 is that the natural frequency and damping of the controller were adjusted 
to allow reasonable computing times per re-entry solution. The numerical 
integration time interval required to produce an accurate solution is a function 
of the latent roots of the equations of motion. Neither aerodynamic character¬ 
istics of the vehicle nor the vehicle's mode of motion (as programmed) results 
in frequencies greater than about 0. 1 radians per second; therefore, the highest 
frequencies would be contributed by the controller dynamics. After this condition 
was recognized, the controller characteristics were changed from the "ideal" 

^ 7 i= 3 rad/sec and !f = 0.7 to a more reasonable &n= 0.3 rad/sec and S’ = 0.7. 
With this lower frequency, the Runge-Kutta integration method produced accurate 
results using a time interval of five seconds. 
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D.2 NUMERICAL INTEGRATION 
D. 2. 1 Introduction 


Twelve time rates of change of both linear and angular position and velocity 
parameters of the vehicle are simultaneously integrated once each computation 
cycle using some numerical method. These parameters are given by Equations 
A. 37, A. 38, A. 48, A.49, A. 50, and A. 55 or A. 55s. The equations are evaluated 
using position and velocity variables assumed given at some time (t ) and thence 
numerically integrated with a time interval or step size (At) to obtain the posi¬ 
tion and velocity variables at the later time t + At . This cycle of operations, 
with the auxiliary computations that may be necessary during each cycle, is 
repeated as long as desired. 

The above comments are applicable to the use of any numerical integration 
method (though specifically referred to the problem of this report). However, 
there is no mention above that a number of numerical integration methods are 
available, each of whose application in a specific situation is different from the 
others. Taylor's Series, the fourth order Runge-Kutta equations, an open equation, 
a closed equation, and an open-closed combination are examples of methods which 
require different mechanics of application although in some cases, the differences 
may be small. Similarly, there is no mention of characteristics of a method 
(per-step error*, propagation of per-step error**, and stable convergence***). 
From the user's point of view, the concept of stable convergence is not of prime 
importance, for if he knows the method is convergent, he need not further concern 
himself with it. Per-step error and its-propagation are, however, very important. 
In fact, these concepts have direct bearing upon whether the integration method may 
be satisfactorily used for real time operations. The studies performed by various 
workers (notably H- J. Gray, Jr.) at the Moore School of Electrical Engineering 
at the University of Pennsylvania were directed at obtaining suitably convergent 
integration formulas with favorable propagation characteristics. A useful result 
of these studies is the so-called stability chart (see Reference D. 3) which provides 
a graphical determination of propagation characteristics. 


* Per-step error as used here includes error from two sources: 

1. Truncation error of the numerical integration formula; that is, the 
formula results in a value which is slightly different from the true 
value. The numerical formulas are generally based upon a finite 
number of terms of a Taylor's or Maclaurin's series; thus the 
infinite series has been "truncated". 

2. Round-off error due to the limited number of digits carried by the 
computer. 

Both errors occur each time the formula is used. 

** Propagation of per-step error concerns how the error obtained at one step 
affects the solution in subsequent steps. 

*** Stable convergence is obtained if, for &t —►O at a point, the value of the 
solution at that point converges to the true solution, providing of course that the 
given data is true. 
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If a method is stably convergent, by the precise definition of Reference 
D.2, the propagation of per-step error is a function of the characteristics of 
the given differential equations and the time interval. The characteristics of 
the equations do, in fact, change with time and, since little control can be had 
over the equations without restricting their solution, the only real freedom 
lies in choosing some optimum time interval. 

Per-step error has been defined previously as the sum of the integra¬ 
tion formula truncation error and the number round-off error. For the Runge- 
Kutta formula discussed in Section D.2.2, the per-step truncation error is 

S T *»*, (*t) r (D.l) 

where At is the integration step-size or time interval and ^ is a constant 
(assumed) which is proportional to the fifth derivative of the solution. The 
per-step round-off error can be represented simply as a constant, , 
which is an average of the possible values of the round-off error and is inde¬ 
pendent of the characteristics of the equations to be solved: 

€ m ~ i x (D. 2) 

Equation D. 1 shows that per-step truncation error, £ r , depends heavily on 
time interval, while Equation D.2 shows that per-step round-off error, £* < 
is independent of time interval. 

Now, over a given length of time history, T , the total accumulated 
per-step error is the sum of accumulated truncation and round-off errors: 




(D.3) 


where T/At is the number of integration steps. The sketch below depicts 
£ , £ m , and £ r for two frequencies, *4 and * ft*, . The one £ m 

curve holds for both frequencies. 


£ 

Per-step error 
accumulated 
over time history 
of length T 
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Ideally, a solution should be obtained by using the time interval resulting in 
minimum accumulated per-step error: this is the time interval at the so-called 
Rademacher point. However, numerical solution of differential equations must 
be accomplished so that requirements of the highest frequency in the equations 
( 60 z in the sketch) are met. As a consequence, the error in the lower fre¬ 
quency ( ) portions of the solution is considerably greater than optimum. 

In addition, the Rademacher point is not fixed; it depends on the char¬ 
acteristics of the differential equations, specifically the latent roots, which 
in general change during the time history. Thus, the equations themselves 
can thwart effort to locate the optimum time interval. 

To obtain solutions of the equations of motion during the present study, 
two numerical integration methods were studied: a fourth order Runge-Kutta 
method and a fourth order open-closed method*. The original plan was to 
modify the open-closed method so that it would automatically change its own 
time interval (to that value which allows the estimated per-step error to meet 
desired limits) thus optimizing computation time with respect to the pre¬ 
selected limits. This method is workable but, because of various complications, 
its development was never completed. As a result, integration was performed 
using the Runge-Kutta method with a fixed-step size (except when changed by 
an independent computer instruction based on considerations other than per- 
step error). The computation time was made acceptable by employing special 
controller characteristics as described in Section 5. 1.6. It became possible, 
therefore, to use a time interval of 5 sec for the time history solutions. 

Certainly dtr = 5 sec is not capable of correctly computing solutions 
if the normal aircraft high frequency dynamics are to be described. The 
solutions presented in Sections V and VI of this report are intended to be 
basically trajectory time histories; nevertheless, conditions have been selected 
so that some information on high frequency dynamics is available. 

D.2.2 Runge-Kutta Method 

Of the many similar formulas called"Runge-Kutta", a particular 
fourth order formula (the most commonly known) was used to obtain the 
solutions presented in Sections V and VI. The formula is given for a single 
first order differential equation: 

***/ s (*, + 2*1+2*, + 4+) (D.4) 

4 , a 

* * t y (** ** ) 

i 9 (tn+J&t, y„ +±i s ) 

* The open-closed method, also called predictor-corrector, provides an 
estimate of per-step error in terms of the predictor solution and the correc¬ 
tor solution. 


where: 


and where: 
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The Runge-Kutta method was chosen as one of two methods to be used with 
the IBM 704 program because of its excellent propagation and truncation 
error properties and because it is always stably convergent and self-starting 
when properly used. Thus it provided reliable solutions with which the modified 
open-closed method could be compared during its development. A drawback 
of the Runge-Kutta method, however, is the necessity for four evaluations 
of the differential equations for each step to obtain the four k's. 

D.2.3 Open-Closed Method 


The open-closed method was selected because it held promise of 
automatically choosing its own optimum time interval,thus providing minimum 
time history computation time. The method-is stably convergent and has 
truncation error of about the same order as Runge-Kutta. The open-closed 
method has generally worse propagation properties, however, and it must 
be started by using some process such as Taylor's series or the Runge-Kutta 
method. The formulas investigated during the study are given in Reference 


OPEN: y*,, - y n +■£ [«y* - ft if,,., f *ry M . M -ty„.^ <d.5> 

CLOSED: HL, * ft (D . 6) 

with the corresponding error formulas: 


OPEN: 


£° 


zn 

7ZO 




(D.7) 


CLOSED: 


£° m 


/4 

720 



(D.8) 


The open formula "predicts" a value , based upon At and , etc., 

which is then used in the closed formula to result in the final "corrected" 
value whose resultant error is £ c . Now £* can be expressed as: 


£ e • 



) (f-e‘) 


(D.9) 


The first factor of Equation D.9 is found to be 14.2 by using Equations D. 7 and 
D.8. Since the true value of the solution is 

^ ° s (D. 10) 
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it follows that 


£°- £ c - y c - y* (D. ii) 

Equation D.9 can evaluated using the different forms of each factor: 


£ 


y e -if 

/*.Z 


(D. 12) 


Equation D. 12 provides an estimate of the per-step error of the open-closed 
method which is dependent only upon the results of the two formulas, and it 
is used with some desired error limits as a criterion for varying time inter¬ 
val. 


The open-closed method requires starting only if the best accuracy of 
the method is wanted. Thus, any time a discontinuity occurs in the equations 
(e.g. , axis switching) the method should be restarted. This reasoning holds 
true for all open, closed, and open-closed formulas but, actually, if the method 
has favorable propagation characteristics and the best accuracy is not required, 
it is not necessary to start or restart. 

The closed formula (Equation D. 6 ) need not be used each computation 
cycle. In fact, in Reference D.4, it is suggested that the closed formula 
be used every tenth step. This would allow a substantial time saving. 

D. 2.4 Conclusion 


The previous sections have presented some of the problems of numerical 
integration and the methods used during the study to overcome these problems. 
However, the techniques discussed here and others that were found during the 
study (see References D.5 and D. 6 ) do not appear to be suitable for simulator 
use. In addition, the author's feeling is that the O 33 mod Gurk formula (see 
Reference D.7) may not be a suitable formula for simulator use with the equations 
of motion developed in Section III, especially if simulated re-entries are long 
enough for severe error accumulation. Unfortunately, O 33 mod Gurk was not 
tried with typical re-entry solutions, so that the "feeling" has little founding 
on fact. 

D. 3 OBTAINING SOLUTIONS OF EQUATIONS OF MOTION 
D.3.1 Initial Condition Calculations for Cases II and III 


The initial values of trajectory parameters for Cases II and III are 
given in Section 5. 2.1 as: 

it = 300,000 ft 
V = 24,000 ft/sec 
7 = -1 deg 
1 ) - +30 deg 
^ = + 120 deg 

Other conditions which pertain are the earth's rotation with angular velocity, 
_TL = 0.004178 deg/sec, the 15* angle of attack to be held from zero time 
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onward, the zero sideslip to be held from zero time onward, zero initial 
angular velocity in inertial space, and zero bank angle. 

Immediately *p = £ = r* = A = 0 =0, and A = 300,000 ft, -i> = 
+ 30 deg. Then because =0 

+1^0 deg 


also because /S = 0 


& » <* + 9 = +14 deg 

Now U , sir , and o'" must be obtained-a process not so trivial as evaluation 
of the preceding nine variables. 

From the geometry of the problem two vector diagrams can be drawn, 
one in a plane horizontal to the local earth's surface and the second perpendicu¬ 
lar to the first and containing the tangential component ( V T ) of the relative 
velocity vector ( V £ ). 



Sketch A 
( not to scale) 



(not to scale) 
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At altitude A - 300,000 ft, and latitude ^ = 30°, the inertial velocity compo¬ 
nent due to earth rotation is: 

v *n s A >Y» * , 

= 0.00007292 (20.86 x 10 b + 0.3 x 10 b ) cos 30 # 
= 1336.2661 ft/sec 

Then, from Sketch A: 

nr = i/ gf , $/„ 30° 

= 668. 13305 ft/sec 


Sketch B shows that 

V* +(V r + £ os 30°)* - u* V ur 2 
which can be written as a quadratic in : 

Vf +(V„ ctt i°+ V mr *'**■*''* 

where ? 

U* + *r* * V*-sv 2 = 575,553,598 (ft/secp 

The physically realizable solution is 

Vf. = 398.501701 ft/sec 
It is now possible to evaluate the angle TP : 


V * s,n ~* f- jt r T> 

y u* + ur 

¥ = 0.9518 deg 

Then finally U. and ur are evaluated: 

U • *■ ur*’ at (&+$>) - 23, 178.450 ft/sec 

tr • y]u* + ur*' „„ (& v p) = 6189.7513 ft/sec 


To check u , ^ , and ur , find the sum of their squares and compare 
this gum with the square of the given V . The 8 digits given for U , ^tr , 
and ur are not an indication of the precision of the values. 

D. 3. 2 Check Solutions and Accuracy 

A number of solutions were obtained early in the study to check vari¬ 
ous portions of the computer program and to demonstrate the accuracy of the 
numerical integration methods over long time histories. For most of these 
special solutions, exact analytical answers could be obtained to compare 
with the numerical solution. 
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The types of check solutions and their purposes are listed here: 

1. A force and moment free solution was obtained to check the rotational 
equations of motion only. This solution was patterned after the special 
integrable force-and moment-free solution on page 162 of Reference D.8. 

2. A circular orbit with - 300,000 ft and the vehicle's attitude fixed 

in inertial space was obtained. This solution was designed to check 
the axis switching mechanism as the vehicle's attitude relative to earth 
changes. It also showed that the Runge-Kutta method could integrate 
orbits fairly well. The orbit did not close on itself but the error was 
only about -6 ft, nearly a half of which could be accounted for by the 
fact that the initial circular orbital velocity could not be programmed 
precisely enough with the 8 digit capacity of the IBM 704. 

3. A highly eccentric (0.36) orbit was computed which is certainly not 

realistic for the earth of radius /? 6 = 20.86 x 10® ft. The altitude at 

apogee is 390,000 ft and at perigee, -10,860,000 ft. The vehicle was 
given angular velocity. This solution was another accuracy check 
with the vehicle undergoing more general motions. Again, the orbital 
accuracy proved good, but not so good as for the circular orbit in (2) 
above. 

4. A constant C D S/j^, re-entry solution was obtained to check use of 
atmospheric models which,until this solution, were unused. 

5. A solution was obtained with conditions set so that an exact second 

order response to an angle of attack step could be computed. This 
was, of course, designed to check the controller formulation and to 
determine the necessary At for reasonable accuracy. The 
solution for CaJ^ - 3 rad/sec and If =0.7 required a time 

interval of 0.5 sec, an entirely unreasonable value from the point of 
view of computation time per solution. 

Check solutions listed here proved very useful in finding programming errors. 

In addition, the last check solution above (5) showed that with At =0.5 sec 
and a 2500 sec time history, the IBM 704, computing at a rate of 40 steps per 
minute, would require over 2 hours for each time history. Two hours was 
entirely too long when some dozens of solutions were required; therefore, a 
method for shortening the computation time per solution had to be devised. 
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